From a31e97c3e0004eaf6c467eeacc831680bf89a684 Mon Sep 17 00:00:00 2001 From: cvs Date: Tue, 12 Oct 1999 17:24:00 +0000 Subject: [PATCH] Fix several bugs that prevented compilation. git-svn-id: https://svn.dealii.org/trunk@1762 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/source/numerics/vectors.cc | 98 ++++++++++++++-------- 1 file changed, 62 insertions(+), 36 deletions(-) diff --git a/deal.II/deal.II/source/numerics/vectors.cc b/deal.II/deal.II/source/numerics/vectors.cc index a12b946c3e..989123d2dd 100644 --- a/deal.II/deal.II/source/numerics/vectors.cc +++ b/deal.II/deal.II/source/numerics/vectors.cc @@ -50,6 +50,11 @@ void VectorTools::interpolate (const DoFHandler &dof, Vector &vec) { const FiniteElement &fe = dof.get_fe(); + + // use #interpolate# function with + // #VectorFunction# param for system + // elements + Assert (fe.n_components == 1, ExcNotUseful()); DoFHandler::active_cell_iterator cell = dof.begin_active(), endc = dof.end(); @@ -82,6 +87,11 @@ void VectorTools::interpolate (const DoFHandler &dof, { const FiniteElement &fe = dof.get_fe(); + // use #interpolate# function with + // #Function# param for non-system + // elements + Assert (fe.n_components == vectorfunction.n_components, ExcNotUseful()); + DoFHandler::active_cell_iterator cell = dof.begin_active(), endc = dof.end(); @@ -100,7 +110,8 @@ void VectorTools::interpolate (const DoFHandler &dof, fe.get_unit_support_points(unit_support_points); // The following works well - // if #dofs_per_cell<=1# as then + // if #dofs_per_x<=1 (x=vertex,line,cell)# + // as then // the multiple support_points // are placed one after another. @@ -273,8 +284,8 @@ void VectorTools::interpolate (const DoFHandler &dof, template void -VectorTools::interpolate(const DoFHandler &high_dof, - const DoFHandler &low_dof, +VectorTools::interpolate(const DoFHandler &high_dof, + const DoFHandler &low_dof, const FullMatrix &transfer, const Vector &high, Vector &low) @@ -284,8 +295,9 @@ VectorTools::interpolate(const DoFHandler &high_dof, DoFHandler::active_cell_iterator h = high_dof.begin_active(); DoFHandler::active_cell_iterator l = low_dof.begin_active(); + const DoFHandler::cell_iterator endh = high_dof.end(); - for(; h != high_dof.end(); ++h, ++l) + for(; h != endh; ++h, ++l) { h->get_dof_values(high, cell_high); transfer.vmult(cell_low, cell_high); @@ -328,7 +340,10 @@ void VectorTools::project (const DoFHandler &dof, Vector &vec, const bool enforce_zero_boundary, const Quadrature &q_boundary, - const bool project_to_boundary_first) { + const bool project_to_boundary_first) +{ + Assert (dof.get_fe().n_components == 1, ExcNotUseful()); + const FiniteElement &fe = dof.get_fe(); // make up boundary values @@ -421,7 +436,10 @@ template void VectorTools::create_right_hand_side (const DoFHandler &dof, const Quadrature &quadrature, const Function &rhs, - Vector &rhs_vector) { + Vector &rhs_vector) +{ + Assert (dof.get_fe().n_components == 1, ExcNotUseful()); + UpdateFlags update_flags = UpdateFlags(update_q_points | update_JxW_values); SparseMatrix dummy; @@ -444,6 +462,7 @@ void VectorTools::create_right_hand_side (const DoFHandler &dof, + #if deal_II_dimension == 1 template <> @@ -457,7 +476,8 @@ VectorTools<1>::interpolate_boundary_values (const DoFHandler<1> &dof, const FiniteElement<1> &fe = dof.get_fe(); Assert (fe.dofs_per_vertex == 1, ExcInvalidFE()); - + Assert (fe.n_components == 1, ExcInvalidFE()); + // check whether boundary values at the // left boundary of the line are requested if (dirichlet_bc.find(0) != dirichlet_bc.end()) @@ -525,13 +545,9 @@ VectorTools::interpolate_boundary_values (const DoFHandler &dof, ExcInvalidBoundaryIndicator()); const FiniteElement &fe = dof.get_fe(); + Assert (fe.dofs_per_vertex == 1, ExcInvalidFE()); + Assert (fe.n_components == 1, ExcInvalidFE()); - // use two face iterators, since we need - // a DoF-iterator for the dof indices, but - // a Tria-iterator for the fe object - DoFHandler::active_face_iterator face = dof.begin_active_face(), - endf = dof.end_face(); - typename FunctionMap::const_iterator function_ptr; // field to store the indices of dofs @@ -539,6 +555,8 @@ VectorTools::interpolate_boundary_values (const DoFHandler &dof, vector > dof_locations (face_dofs.size(), Point()); vector dof_values (fe.dofs_per_face); + DoFHandler::active_face_iterator face = dof.begin_active_face(), + endf = dof.end_face(); for (; face!=endf; ++face) if ((function_ptr = dirichlet_bc.find(face->boundary_indicator())) != dirichlet_bc.end()) @@ -570,20 +588,19 @@ VectorTools::interpolate_boundary_values (const DoFHandler &dof, ExcInvalidBoundaryIndicator()); const FiniteElement &fe = dof.get_fe(); + Assert (fe.n_components == dirichlet_bc.begin()->second->n_components, + ExcInvalidFE()); - // use two face iterators, since we need - // a DoF-iterator for the dof indices, but - // a Tria-iterator for the fe object - DoFHandler::active_face_iterator face = dof.begin_active_face(), - endf = dof.end_face(); - typename VectorFunctionMap::const_iterator function_ptr; // field to store the indices of dofs - vector face_dofs (fe.dofs_per_face); + vector face_dofs (fe.dofs_per_face, -1); vector > dof_locations (face_dofs.size(), Point()); - vector< Vector > dof_values (fe.dofs_per_face, Vector(fe.n_components)); + vector< Vector > dof_values (fe.dofs_per_face, + Vector(fe.n_components)); + DoFHandler::active_face_iterator face = dof.begin_active_face(), + endf = dof.end_face(); for (; face!=endf; ++face) if ((function_ptr = dirichlet_bc.find(face->boundary_indicator())) != dirichlet_bc.end()) @@ -600,13 +617,9 @@ VectorTools::interpolate_boundary_values (const DoFHandler &dof, // enter into list for (unsigned int i=0; i - index = fe.face_system_to_component_index(i); - const double s = dof_values[i](index.first); - boundary_values[face_dofs[i]] = s; - } - } + boundary_values[face_dofs[i]] + = dof_values[i](fe.face_system_to_component_index(i).first); + }; } @@ -617,6 +630,8 @@ VectorTools::project_boundary_values (const DoFHandler &dof, const FunctionMap &boundary_functions, const Quadrature &q, map &boundary_values) { + Assert (dof.get_fe().n_components == 1, ExcInvalidFE()); + vector dof_to_boundary_mapping; dof.map_dof_to_boundary_indices (boundary_functions, dof_to_boundary_mapping); @@ -1104,23 +1119,34 @@ VectorTools::integrate_difference (const DoFHandler &dof, }; }; + + template void -VectorTools::subtract_mean_value(Vector& v, const - bit_vector& p_select) +VectorTools::subtract_mean_value(Vector &v, + const vector &p_select) { unsigned int n = v.size(); Assert(n == p_select.size(), ExcDimensionMismatch(n, p_select.size())); - double s = 0; - for (unsigned int i=0;i; -- 2.39.5