]> https://gitweb.dealii.org/ - dealii.git/commitdiff
step-50: better matrix construction
authorTimo Heister <timo.heister@gmail.com>
Sat, 25 May 2013 15:24:31 +0000 (15:24 +0000)
committerTimo Heister <timo.heister@gmail.com>
Sat, 25 May 2013 15:24:31 +0000 (15:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@29597 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-50/step-50.cc

index d446b0a2469e62497eed0cb447f43c24cb442bd7..abca04f105575435e15b6b8271a1d2e78fe5db41 100644 (file)
@@ -135,7 +135,6 @@ namespace Step50
     typedef TrilinosWrappers::SparseMatrix matrix_t;
     typedef TrilinosWrappers::MPI::Vector vector_t;
 
-    SparsityPattern      sparsity_pattern;
     matrix_t system_matrix;
 
     // We need an additional object for the
@@ -191,7 +190,6 @@ namespace Step50
     // transpose of the other, and so we only
     // have to build one; we choose the one
     // from coarse to fine.
-    MGLevelObject<SparsityPattern>       mg_sparsity_patterns;
     MGLevelObject<matrix_t> mg_matrices;
     MGLevelObject<matrix_t> mg_interface_matrices;
     MGConstrainedDoFs                    mg_constrained_dofs;
@@ -331,11 +329,6 @@ namespace Step50
               << mg_dof_handler.n_dofs(l);
     deallog  << std::endl;
 
-    sparsity_pattern.reinit (mg_dof_handler.n_dofs(),
-                             mg_dof_handler.n_dofs(),
-                             mg_dof_handler.max_couplings_between_dofs());
-    DoFTools::make_sparsity_pattern (mg_dof_handler, sparsity_pattern);
-
 
     //solution.reinit (mg_dof_handler.n_dofs());
     //system_rhs.reinit (mg_dof_handler.n_dofs());
@@ -375,9 +368,10 @@ namespace Step50
                                               constraints);
     constraints.close ();
     hanging_node_constraints.close ();
-    constraints.condense (sparsity_pattern);
-    sparsity_pattern.compress();
-    system_matrix.reinit (sparsity_pattern);
+
+    CompressedSimpleSparsityPattern csp(mg_dof_handler.n_dofs(), mg_dof_handler.n_dofs());
+    DoFTools::make_sparsity_pattern (mg_dof_handler, csp, constraints);
+    system_matrix.reinit (mg_dof_handler.locally_owned_dofs(), csp, MPI_COMM_WORLD, true);
 
     // The multigrid constraints have to be
     // initialized. They need to know about
@@ -409,7 +403,6 @@ namespace Step50
     mg_interface_matrices.clear ();
     mg_matrices.resize(0, n_levels-1);
     mg_matrices.clear ();
-    mg_sparsity_patterns.resize(0, n_levels-1);
 
     // Now, we have to provide a matrix on each
     // level. To this end, we first use the
@@ -444,10 +437,10 @@ namespace Step50
                    mg_dof_handler.n_dofs(level));
         MGTools::make_sparsity_pattern(mg_dof_handler, csp, level);
 
-        mg_sparsity_patterns[level].copy_from (csp);
-
-        mg_matrices[level].reinit(mg_sparsity_patterns[level]);
-        mg_interface_matrices[level].reinit(mg_sparsity_patterns[level]);
+        mg_matrices[level].reinit(mg_dof_handler.locally_owned_mg_dofs(level),
+            mg_dof_handler.locally_owned_mg_dofs(level),
+            csp,
+            MPI_COMM_WORLD, true);
       }
   }
 

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.