]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Direct-initialize matrix from dsp.
authorWolfgang Bangerth <bangerth@colostate.edu>
Sun, 2 Jul 2023 14:36:15 +0000 (08:36 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Sun, 2 Jul 2023 14:36:15 +0000 (08:36 -0600)
examples/step-86/step-86.cc

index bad3e2f8792d930f13d2be81269dfd885fc81b14..6981ae53d48ca6c2f7c939040602064b42f79358 100644 (file)
@@ -28,7 +28,6 @@
 #include <deal.II/lac/petsc_vector.h>
 #include <deal.II/lac/petsc_sparse_matrix.h>
 #include <deal.II/lac/petsc_solver.h>
-#include <deal.II/lac/sparsity_pattern.h>
 #include <deal.II/lac/petsc_precondition.h>
 #include <deal.II/lac/affine_constraints.h>
 #include <deal.II/grid/tria.h>
@@ -74,7 +73,6 @@ namespace Step86
 
     AffineConstraints<double> constraints;
 
-    SparsityPattern                  sparsity_pattern;
     PETScWrappers::MPI::SparseMatrix mass_matrix;
     PETScWrappers::MPI::SparseMatrix laplace_matrix;
     PETScWrappers::MPI::SparseMatrix system_matrix;
@@ -195,17 +193,11 @@ namespace Step86
                                     dsp,
                                     constraints,
                                     /*keep_constrained_dofs = */ true);
-    sparsity_pattern.copy_from(dsp);
-
-    mass_matrix.reinit(dof_handler.locally_owned_dofs(),
-                       sparsity_pattern,
-                       MPI_COMM_SELF);
-    laplace_matrix.reinit(dof_handler.locally_owned_dofs(),
-                          sparsity_pattern,
-                          MPI_COMM_SELF);
-    system_matrix.reinit(dof_handler.locally_owned_dofs(),
-                         sparsity_pattern,
-                         MPI_COMM_SELF);
+
+    // directly initialize from dsp, no need for the regular sparsity pattern:
+    mass_matrix.reinit(dof_handler.locally_owned_dofs(), dsp, MPI_COMM_SELF);
+    laplace_matrix.reinit(dof_handler.locally_owned_dofs(), dsp, MPI_COMM_SELF);
+    system_matrix.reinit(dof_handler.locally_owned_dofs(), dsp, MPI_COMM_SELF);
 
     MatrixCreator::create_mass_matrix(dof_handler,
                                       QGauss<dim>(fe.degree + 1),

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.