]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Simplify initialization of block objects in the tutorial. 13109/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 23 Dec 2021 02:11:55 +0000 (19:11 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 23 Dec 2021 02:11:55 +0000 (19:11 -0700)
examples/step-22/step-22.cc
examples/step-31/step-31.cc
examples/step-43/step-43.cc
examples/step-44/step-44.cc
examples/step-79/step-79.cc
source/lac/trilinos_block_vector.cc

index 1ab4c866c77fa1545fbb4a0223b3815624065b30..c24fa9e0968125da42aec6191a82dd3aa2b6e7d2 100644 (file)
@@ -535,17 +535,9 @@ namespace Step22
     // <code>dsp</code> will be released once the information has been copied to
     // <code>sparsity_pattern</code>.
     {
-      BlockDynamicSparsityPattern dsp(2, 2);
-
-      dsp.block(0, 0).reinit(n_u, n_u);
-      dsp.block(1, 0).reinit(n_p, n_u);
-      dsp.block(0, 1).reinit(n_u, n_p);
-      dsp.block(1, 1).reinit(n_p, n_p);
-
-      dsp.collect_sizes();
+      BlockDynamicSparsityPattern dsp(dofs_per_block, dofs_per_block);
 
       Table<2, DoFTools::Coupling> coupling(dim + 1, dim + 1);
-
       for (unsigned int c = 0; c < dim + 1; ++c)
         for (unsigned int d = 0; d < dim + 1; ++d)
           if (!((c == dim) && (d == dim)))
@@ -560,17 +552,10 @@ namespace Step22
     }
 
     {
-      BlockDynamicSparsityPattern preconditioner_dsp(2, 2);
-
-      preconditioner_dsp.block(0, 0).reinit(n_u, n_u);
-      preconditioner_dsp.block(1, 0).reinit(n_p, n_u);
-      preconditioner_dsp.block(0, 1).reinit(n_u, n_p);
-      preconditioner_dsp.block(1, 1).reinit(n_p, n_p);
-
-      preconditioner_dsp.collect_sizes();
+      BlockDynamicSparsityPattern preconditioner_dsp(dofs_per_block,
+                                                     dofs_per_block);
 
       Table<2, DoFTools::Coupling> preconditioner_coupling(dim + 1, dim + 1);
-
       for (unsigned int c = 0; c < dim + 1; ++c)
         for (unsigned int d = 0; d < dim + 1; ++d)
           if (((c == dim) && (d == dim)))
@@ -593,15 +578,8 @@ namespace Step22
     system_matrix.reinit(sparsity_pattern);
     preconditioner_matrix.reinit(preconditioner_sparsity_pattern);
 
-    solution.reinit(2);
-    solution.block(0).reinit(n_u);
-    solution.block(1).reinit(n_p);
-    solution.collect_sizes();
-
-    system_rhs.reinit(2);
-    system_rhs.block(0).reinit(n_u);
-    system_rhs.block(1).reinit(n_p);
-    system_rhs.collect_sizes();
+    solution.reinit(dofs_per_block);
+    system_rhs.reinit(dofs_per_block);
   }
 
 
index 7d6463e33c69430819bf8194f9b9f951ef714462..380e9185191c28c24cf12dca6306d7abc699c482 100644 (file)
@@ -870,6 +870,7 @@ namespace Step31
     const unsigned int n_u = stokes_dofs_per_block[0],
                        n_p = stokes_dofs_per_block[1],
                        n_T = temperature_dof_handler.n_dofs();
+    const std::vector<unsigned int> stokes_block_sizes = {n_u, n_p};
 
     std::cout << "Number of active cells: " << triangulation.n_active_cells()
               << " (on " << triangulation.n_levels() << " levels)" << std::endl
@@ -919,14 +920,7 @@ namespace Step31
     {
       stokes_matrix.clear();
 
-      BlockDynamicSparsityPattern dsp(2, 2);
-
-      dsp.block(0, 0).reinit(n_u, n_u);
-      dsp.block(0, 1).reinit(n_u, n_p);
-      dsp.block(1, 0).reinit(n_p, n_u);
-      dsp.block(1, 1).reinit(n_p, n_p);
-
-      dsp.collect_sizes();
+      BlockDynamicSparsityPattern dsp(stokes_block_sizes, stokes_block_sizes);
 
       Table<2, DoFTools::Coupling> coupling(dim + 1, dim + 1);
 
@@ -948,14 +942,7 @@ namespace Step31
       Mp_preconditioner.reset();
       stokes_preconditioner_matrix.clear();
 
-      BlockDynamicSparsityPattern dsp(2, 2);
-
-      dsp.block(0, 0).reinit(n_u, n_u);
-      dsp.block(0, 1).reinit(n_u, n_p);
-      dsp.block(1, 0).reinit(n_p, n_u);
-      dsp.block(1, 1).reinit(n_p, n_p);
-
-      dsp.collect_sizes();
+      BlockDynamicSparsityPattern dsp(stokes_block_sizes, stokes_block_sizes);
 
       Table<2, DoFTools::Coupling> coupling(dim + 1, dim + 1);
       for (unsigned int c = 0; c < dim + 1; ++c)
index 7af90f333e98c9fff2fef63b252a166e6692131b..e29f03443c3ffb9300c1a3614f5a9199fc308d49 100644 (file)
@@ -695,6 +695,7 @@ namespace Step43
     const unsigned int n_u = darcy_dofs_per_block[0],
                        n_p = darcy_dofs_per_block[1],
                        n_s = saturation_dof_handler.n_dofs();
+    const std::vector<unsigned int> darcy_block_sizes = {n_u, n_p};
 
     std::cout << "Number of active cells: " << triangulation.n_active_cells()
               << " (on " << triangulation.n_levels() << " levels)" << std::endl
@@ -705,17 +706,9 @@ namespace Step43
     {
       darcy_matrix.clear();
 
-      BlockDynamicSparsityPattern dsp(2, 2);
-
-      dsp.block(0, 0).reinit(n_u, n_u);
-      dsp.block(0, 1).reinit(n_u, n_p);
-      dsp.block(1, 0).reinit(n_p, n_u);
-      dsp.block(1, 1).reinit(n_p, n_p);
-
-      dsp.collect_sizes();
+      BlockDynamicSparsityPattern dsp(darcy_block_sizes, darcy_block_sizes);
 
       Table<2, DoFTools::Coupling> coupling(dim + 1, dim + 1);
-
       for (unsigned int c = 0; c < dim + 1; ++c)
         for (unsigned int d = 0; d < dim + 1; ++d)
           if (!((c == dim) && (d == dim)))
@@ -735,14 +728,7 @@ namespace Step43
       Mp_preconditioner.reset();
       darcy_preconditioner_matrix.clear();
 
-      BlockDynamicSparsityPattern dsp(2, 2);
-
-      dsp.block(0, 0).reinit(n_u, n_u);
-      dsp.block(0, 1).reinit(n_u, n_p);
-      dsp.block(1, 0).reinit(n_p, n_u);
-      dsp.block(1, 1).reinit(n_p, n_p);
-
-      dsp.collect_sizes();
+      BlockDynamicSparsityPattern dsp(darcy_block_sizes, darcy_block_sizes);
 
       Table<2, DoFTools::Coupling> coupling(dim + 1, dim + 1);
       for (unsigned int c = 0; c < dim + 1; ++c)
@@ -773,23 +759,19 @@ namespace Step43
       saturation_matrix.reinit(dsp);
     }
 
-    std::vector<IndexSet> darcy_partitioning(2);
-    darcy_partitioning[0] = complete_index_set(n_u);
-    darcy_partitioning[1] = complete_index_set(n_p);
+    const std::vector<IndexSet> darcy_partitioning = {complete_index_set(n_u),
+                                                      complete_index_set(n_p)};
+
     darcy_solution.reinit(darcy_partitioning, MPI_COMM_WORLD);
-    darcy_solution.collect_sizes();
 
     last_computed_darcy_solution.reinit(darcy_partitioning, MPI_COMM_WORLD);
-    last_computed_darcy_solution.collect_sizes();
 
     second_last_computed_darcy_solution.reinit(darcy_partitioning,
                                                MPI_COMM_WORLD);
-    second_last_computed_darcy_solution.collect_sizes();
 
     darcy_rhs.reinit(darcy_partitioning, MPI_COMM_WORLD);
-    darcy_rhs.collect_sizes();
 
-    IndexSet saturation_partitioning = complete_index_set(n_s);
+    const IndexSet saturation_partitioning = complete_index_set(n_s);
     saturation_solution.reinit(saturation_partitioning, MPI_COMM_WORLD);
     old_saturation_solution.reinit(saturation_partitioning, MPI_COMM_WORLD);
     old_old_saturation_solution.reinit(saturation_partitioning, MPI_COMM_WORLD);
index 0a5e7c9e4aadd556b5054e1456c45b6174acd4fb..082ddf6fedbf0c4154dd31401183d7f6bb25c18e 100644 (file)
@@ -1449,24 +1449,7 @@ namespace Step44
     // Setup the sparsity pattern and tangent matrix
     tangent_matrix.clear();
     {
-      const types::global_dof_index n_dofs_u = dofs_per_block[u_dof];
-      const types::global_dof_index n_dofs_p = dofs_per_block[p_dof];
-      const types::global_dof_index n_dofs_J = dofs_per_block[J_dof];
-
-      BlockDynamicSparsityPattern dsp(n_blocks, n_blocks);
-
-      dsp.block(u_dof, u_dof).reinit(n_dofs_u, n_dofs_u);
-      dsp.block(u_dof, p_dof).reinit(n_dofs_u, n_dofs_p);
-      dsp.block(u_dof, J_dof).reinit(n_dofs_u, n_dofs_J);
-
-      dsp.block(p_dof, u_dof).reinit(n_dofs_p, n_dofs_u);
-      dsp.block(p_dof, p_dof).reinit(n_dofs_p, n_dofs_p);
-      dsp.block(p_dof, J_dof).reinit(n_dofs_p, n_dofs_J);
-
-      dsp.block(J_dof, u_dof).reinit(n_dofs_J, n_dofs_u);
-      dsp.block(J_dof, p_dof).reinit(n_dofs_J, n_dofs_p);
-      dsp.block(J_dof, J_dof).reinit(n_dofs_J, n_dofs_J);
-      dsp.collect_sizes();
+      BlockDynamicSparsityPattern dsp(dofs_per_block, dofs_per_block);
 
       // The global system matrix initially has the following structure
       // @f{align*}
@@ -1511,10 +1494,7 @@ namespace Step44
 
     // We then set up storage vectors
     system_rhs.reinit(dofs_per_block);
-    system_rhs.collect_sizes();
-
     solution_n.reinit(dofs_per_block);
-    solution_n.collect_sizes();
 
     // ...and finally set up the quadrature
     // point history:
index 8018bd81cf00680f77ca97e9333db1d8a725f21a..e1dd2e4dc57e57a16e924222de7a71d792bb9c06 100644 (file)
@@ -441,11 +441,7 @@ namespace SAND
     const std::vector<BlockVector<double>::size_type> block_sizes = {
       n_p, n_u, n_p, n_u, n_p, n_p, n_p, n_p, n_p};
 
-    BlockDynamicSparsityPattern dsp(9, 9);
-    for (unsigned int k = 0; k < 9; ++k)
-      for (unsigned int j = 0; j < 9; ++j)
-        dsp.block(j, k).reinit(block_sizes[j], block_sizes[k]);
-    dsp.collect_sizes();
+    BlockDynamicSparsityPattern dsp(block_sizes, block_sizes);
 
 
     // The bulk of the function is in setting up which of these
index 628e5da397f9e64b098dd29aa559b80b74418b34..25b25c1d621fbd654711f0912c5df97ef1eae735 100644 (file)
@@ -91,6 +91,8 @@ namespace TrilinosWrappers
       collect_sizes();
     }
 
+
+
     void
     BlockVector::reinit(const std::vector<IndexSet> &parallel_partitioning,
                         const std::vector<IndexSet> &ghost_values,
@@ -119,6 +121,7 @@ namespace TrilinosWrappers
     }
 
 
+
     void
     BlockVector::reinit(const BlockVector &v, const bool omit_zeroing_entries)
     {

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.