]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Minor fixes
authorAlexander Grayver <agrayver@erdw.ethz.ch>
Fri, 9 Oct 2015 16:07:20 +0000 (18:07 +0200)
committerAlexander Grayver <agrayver@erdw.ethz.ch>
Fri, 9 Oct 2015 16:48:54 +0000 (18:48 +0200)
source/dofs/dof_tools_constraints.cc
tests/distributed_grids/intergrid_transfer_representation_parallel.cc
tests/distributed_grids/intergrid_transfer_representation_parallel.mpirun=3.output [moved from tests/distributed_grids/intergrid_transfer_representation_parallel.output with 100% similarity]

index 81128949c821633c6386d79facea431661fa4421..cf74c4b2c6674e1fa70a90c8415f051d713c53db 100644 (file)
@@ -2623,7 +2623,8 @@ namespace DoFTools
                                     const unsigned int                                      coarse_component,
                                     const FiniteElement<dim,spacedim>                      &coarse_fe,
                                     const std::vector<types::global_dof_index>             &weight_mapping,
-                                    std::vector<std::map<types::global_dof_index, float> > &weights)
+                                    std::vector<std::map<types::global_dof_index, float> > &weights,
+                                    bool is_called_in_parallel)
       {
         unsigned int pos = 0;
         for (unsigned int local_dof=0; local_dof<copy_data.dofs_per_cell; ++local_dof)
@@ -2667,9 +2668,12 @@ namespace DoFTools
                         weights[wi][wj] = copy_data.global_parameter_representation[pos](i);
                       }
                   }
-              //else
-              //  Assert (copy_data.global_parameter_representation[pos](i) == 0,
-              //          ExcInternalError());
+                else
+                  // Note that when this function operates with distributed
+                  // fine grid, this assertion is switched off since the
+                  // condition does not necessarily hold
+                  Assert (copy_data.global_parameter_representation[pos](i) == 0 || is_called_in_parallel,
+                          ExcInternalError());
               ++pos;
             }
 
@@ -2704,6 +2708,7 @@ namespace DoFTools
 
         copy_data.global_parameter_representation.resize(n_interesting_dofs);
 
+        bool is_called_in_parallel = false;
         for (size_t i=0; i<copy_data.global_parameter_representation.size(); ++i)
           {
 #ifdef DEAL_II_WITH_MPI
@@ -2714,6 +2719,7 @@ namespace DoFTools
                   dynamic_cast<const typename dealii::parallel::Triangulation<dim, spacedim>&>
                   (coarse_to_fine_grid_map.get_destination_grid().get_tria ());
                 communicator = tria.get_communicator ();
+                is_called_in_parallel = true;
               }
             catch (std::bad_cast &exp)
               {
@@ -2750,7 +2756,8 @@ namespace DoFTools
                                         coarse_component,
                                         std_cxx11::cref(coarse_grid.get_fe()),
                                         std_cxx11::cref(weight_mapping),
-                                        std_cxx11::ref(weights)),
+                                        std_cxx11::ref(weights),
+                                        is_called_in_parallel),
                         scratch,
                         copy_data);
       }
index 3b4a7c8c3811bb790443de69b1b006b0e7172bcf..03c63ed11e0ed5591c857ea80e0994e45cd27d2e 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------\r
 //\r
-// Copyright (C) 2013 - 2014 by the deal.II authors\r
+// Copyright (C) 2013 - 2015 by the deal.II authors\r
 //\r
 // This file is part of the deal.II library.\r
 //\r

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.