]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add support for numbering of boundary (for integration on the boundary).
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 22 May 1998 07:39:16 +0000 (07:39 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 22 May 1998 07:39:16 +0000 (07:39 +0000)
git-svn-id: https://svn.dealii.org/trunk@330 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/source/dofs/dof_handler.cc

index fbb3c795dcf94555a77a3ee85daf22ba889885b7..dc5049ac12d795216446d5c43397138a159ef20e 100644 (file)
@@ -725,8 +725,8 @@ void DoFHandler<dim>::distribute_dofs (const FiniteElementBase<dim> &fe) {
 
 
 
-int DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator   &cell,
-                                           unsigned int            next_free_dof) {
+unsigned int DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator &cell,
+                                                    unsigned int          next_free_dof) {
 
                                   // distribute dofs of vertices
   for (unsigned int v=0; v<GeometryInfo<1>::vertices_per_cell; ++v)
@@ -776,8 +776,8 @@ int DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator   &cell,
 
 
 
-int DoFHandler<2>::distribute_dofs_on_cell (active_cell_iterator   &cell,
-                                           unsigned int            next_free_dof) {
+unsigned int DoFHandler<2>::distribute_dofs_on_cell (active_cell_iterator &cell,
+                                                    unsigned int          next_free_dof) {
   if (selected_fe->dofs_per_vertex > 0)
                                     // number dofs on vertices
     for (unsigned int vertex=0; vertex<GeometryInfo<2>::vertices_per_cell; ++vertex)
@@ -824,7 +824,7 @@ int DoFHandler<2>::distribute_dofs_on_cell (active_cell_iterator   &cell,
 
 template <int dim>
 void DoFHandler<dim>::renumber_dofs (const RenumberingMethod method,
-                                    bool use_constraints,
+                                    const bool use_constraints,
                                     const vector<int> &starting_points) {
                                   // make the connection graph
   dSMatrixStruct sparsity (n_dofs(), max_couplings_between_dofs());
@@ -1185,7 +1185,7 @@ void DoFHandler<dim>::make_boundary_sparsity_pattern (const vector<int> &dof_to_
        
                                         // make sparsity pattern for this cell
        for (unsigned int i=0; i<total_dofs; ++i)
-         for (unsigned int j=0; j<total_dofs; ++j)
+         for (unsigned int j=0; j<total_dofs; ++j) 
            sparsity.add (dof_to_boundary_mapping[dofs_on_this_face[i]],
                          dof_to_boundary_mapping[dofs_on_this_face[j]]);
       };
@@ -1625,7 +1625,7 @@ void DoFHandler<dim>::map_dof_to_boundary_indices (vector<int> &mapping) const {
 
   mapping.clear ();
   mapping.insert (mapping.end(), n_dofs(), -1);
-
+  
   const unsigned int dofs_per_face = selected_fe->dofs_per_face;
   vector<int> dofs_on_face(dofs_per_face);
   int next_boundary_index = 0;
@@ -1665,6 +1665,10 @@ void DoFHandler<dim>::map_dof_to_boundary_indices (const FunctionMap &boundary_i
   mapping.clear ();
   mapping.insert (mapping.end(), n_dofs(), -1);
 
+                                  // return if there is nothing to do
+  if (boundary_indicators.size() == 0)
+    return;
+  
   const unsigned int dofs_per_face = selected_fe->dofs_per_face;
   vector<int> dofs_on_face(dofs_per_face);
   int next_boundary_index = 0;

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.