]> https://gitweb.dealii.org/ - dealii.git/commitdiff
DoFTools::extract_constant_modes() and extract_rigid_body_modes(): make ComponentMask... 17637/head
authorPeter Munch <peterrmuench@gmail.com>
Fri, 30 Aug 2024 15:43:36 +0000 (17:43 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Sun, 1 Sep 2024 20:40:07 +0000 (22:40 +0200)
examples/step-42/step-42.cc
examples/step-90/step-90.cc
include/deal.II/dofs/dof_tools.h
source/dofs/dof_tools.cc
tests/mpi/extract_rigid_body_modes_01.cc

index c4ef650a1a4a50a66c41dc0c0f9d41a23b7269ef..f601e1baf89e99c1c5b73f1b8d77e9e892475452 100644 (file)
@@ -1623,7 +1623,7 @@ namespace Step42
       TimerOutput::Scope t(computing_timer, "Solve: setup preconditioner");
 
       const std::vector<std::vector<bool>> constant_modes =
-        DoFTools::extract_constant_modes(dof_handler, ComponentMask());
+        DoFTools::extract_constant_modes(dof_handler);
 
       TrilinosWrappers::PreconditionAMG::AdditionalData additional_data;
       additional_data.constant_modes        = constant_modes;
index c14ae9b33199d2b333a3702221526e1347515bbf..d09a3b0a57572650ae4b9c6fd3ee1a0c294d4f3b 100644 (file)
@@ -947,7 +947,7 @@ namespace Step90
         const unsigned int max_iterations = 500;
         SolverControl      solver_control(max_iterations, relative_error);
         const std::vector<std::vector<bool>> constant_modes =
-          DoFTools::extract_constant_modes(dof_handler, ComponentMask());
+          DoFTools::extract_constant_modes(dof_handler);
         TrilinosWrappers::PreconditionAMG preconditioner_stiffness;
         TrilinosWrappers::PreconditionAMG::AdditionalData Amg_data;
         Amg_data.constant_modes        = constant_modes;
index 0c25937c1d91046e81666633d2e624b116dbb15f..422136e0ad31809d99e67d275df1fdd8036a427d 100644 (file)
@@ -1373,7 +1373,7 @@ namespace DoFTools
   template <int dim, int spacedim>
   std::vector<std::vector<bool>>
   extract_constant_modes(const DoFHandler<dim, spacedim> &dof_handler,
-                         const ComponentMask             &component_mask);
+                         const ComponentMask             &component_mask = {});
 
   /**
    * Same as above.
@@ -1394,7 +1394,7 @@ namespace DoFTools
   std::vector<std::vector<bool>>
   extract_level_constant_modes(const unsigned int               level,
                                const DoFHandler<dim, spacedim> &dof_handler,
-                               const ComponentMask             &component_mask);
+                               const ComponentMask &component_mask = {});
 
   /**
    * Same as above.
@@ -1419,7 +1419,7 @@ namespace DoFTools
   std::vector<std::vector<double>>
   extract_rigid_body_modes(const Mapping<dim, spacedim>    &mapping,
                            const DoFHandler<dim, spacedim> &dof_handler,
-                           const ComponentMask             &component_mask);
+                           const ComponentMask &component_mask = {});
 
   /**
    * Same as above but for multigrid levels.
@@ -1429,7 +1429,7 @@ namespace DoFTools
   extract_level_rigid_body_modes(const unsigned int               level,
                                  const Mapping<dim, spacedim>    &mapping,
                                  const DoFHandler<dim, spacedim> &dof_handler,
-                                 const ComponentMask &component_mask);
+                                 const ComponentMask &component_mask = {});
   /** @} */
 
   /**
index 497dc3bcc8dbf0d598eb1faf2853bcd26f0c2128..4178e691593b1a811bd1f0ba437997ded3f73901 100644 (file)
@@ -1504,7 +1504,6 @@ namespace DoFTools
                              const unsigned int               mg_level)
     {
       AssertDimension(dim, spacedim);
-      AssertDimension(component_mask.n_selected_components(), dim);
 
       constexpr unsigned int n_modes = RigidBodyMotion<dim>::n_modes;
 
index 78e247c5cb3e78ea684d4f805f7e584309337670..e9288a19f01be2cd9645ac674e5e698f8924a917 100644 (file)
@@ -57,10 +57,9 @@ test()
   deallog << "Total dofs=" << dofh.n_dofs() << std::endl;
 
   // extract constant modes and print
-  ComponentMask mask(fe.n_components(), true);
 
   std::vector<std::vector<double>> constant_modes =
-    DoFTools::extract_rigid_body_modes(mapping, dofh, mask);
+    DoFTools::extract_rigid_body_modes(mapping, dofh);
 
   for (unsigned int i = 0; i < constant_modes.size(); ++i)
     {

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.