From: wolf Date: Sun, 29 Apr 2001 16:01:11 +0000 (+0000) Subject: Improve algorithm in create_boundary_mass_matrix. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f9a441601cb9ae7b194acace720de3fc53729c78;p=dealii-svn.git Improve algorithm in create_boundary_mass_matrix. git-svn-id: https://svn.dealii.org/trunk@4499 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/numerics/matrices.cc b/deal.II/deal.II/source/numerics/matrices.cc index 2c1f2ca4e3..0f4d6dfa40 100644 --- a/deal.II/deal.II/source/numerics/matrices.cc +++ b/deal.II/deal.II/source/numerics/matrices.cc @@ -312,19 +312,6 @@ void MatrixCreator::create_mass_matrix (const DoFHandler &dof, #if deal_II_dimension == 1 -template <> -void MatrixCreator<1>::create_boundary_mass_matrix (const DoFHandler<1> &, - const Quadrature<0> &, - SparseMatrix &, - const FunctionMap<1>::type&, - Vector &, - std::vector &, - const Function<1> *) -{ - Assert (false, ExcNotImplemented()); -}; - - template <> void MatrixCreator<1>::create_boundary_mass_matrix (const Mapping<1> &, const DoFHandler<1> &, @@ -335,6 +322,9 @@ void MatrixCreator<1>::create_boundary_mass_matrix (const Mapping<1> &, std::vector &, const Function<1> *) { + // what would that be in 1d? the + // identity matrix on tehe boundary + // dofs? Assert (false, ExcNotImplemented()); }; @@ -459,30 +449,34 @@ create_boundary_mass_matrix_1 (const Mapping &mapping, // the name std::vector coefficient_values_scalar (fe_values.n_quadrature_points); std::vector > coefficient_values_system (fe_values.n_quadrature_points, - Vector(n_components)); + Vector(n_components)); std::vector rhs_values_scalar (fe_values.n_quadrature_points); std::vector > rhs_values_system (fe_values.n_quadrature_points, - Vector(n_components)); + Vector(n_components)); std::vector dofs (dofs_per_cell); std::vector dofs_on_face_vector (dofs_per_face); - std::set dofs_on_face; - DoFHandler::active_cell_iterator cell = range.first; + // for each dof on the cell, have a + // flag whether it is on the face + std::vector dof_is_on_face(dofs_per_cell); + + typename DoFHandler::active_cell_iterator cell = range.first; for (; cell!=range.second; ++cell) for (unsigned int face=0; face::faces_per_cell; ++face) // check if this face is on that part of // the boundary we are interested in - if (boundary_functions.find(cell->face(face)->boundary_indicator()) != boundary_functions.end()) + if (boundary_functions.find(cell->face(face)->boundary_indicator()) != + boundary_functions.end()) { cell_matrix.clear (); cell_vector.clear (); fe_values.reinit (cell, face); - const FullMatrix &values = fe_values.get_shape_values (); - const std::vector &weights = fe_values.get_JxW_values (); + const FullMatrix &values = fe_values.get_shape_values (); + const std::vector &weights = fe_values.get_JxW_values (); if (fe_is_system) // FE has several components @@ -623,9 +617,16 @@ create_boundary_mass_matrix_1 (const Mapping &mapping, cell->get_dof_indices (dofs); cell->face(face)->get_dof_indices (dofs_on_face_vector); - dofs_on_face.clear (); - dofs_on_face.insert (dofs_on_face_vector.begin(), - dofs_on_face_vector.end()); + // check for each of the + // dofs on this cell + // whether it is on the + // face + for (unsigned int i=0; i &mapping, if (fabs(cell_matrix(i,i)) > max_diag_entry) max_diag_entry = fabs(cell_matrix(i,i)); #endif -//TODO: [WB] check this place for more efficient alternatives -//in the innermost loop, we traverse a set twice in "find", but -//this could be made much faster, e.g. by first building a vector -//that already stores the result of find beforehand -// (then remove dofs_on_face altogether) + mutex.acquire (); for (unsigned int i=0; i &mapping, }; for (unsigned int j=0; j