]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
work on step-32 with PETSc (lots missing)
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 Apr 2013 00:08:46 +0000 (00:08 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 Apr 2013 00:08:46 +0000 (00:08 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_unify_linear_algebra@29271 0785d39b-7218-0410-832d-ea1e28bc413d

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

index dd6625b7522fe72275d37684381d3a87e68f97f3..e6c9889a8ef2418a4296172efbd260b5276a5c4d 100644 (file)
@@ -63,7 +63,7 @@
 
 #include <deal.II/lac/abstract_linear_algebra.h>
 
-//#define USE_PETSC_LA
+#define USE_PETSC_LA
 
 namespace LA
 {
@@ -1713,9 +1713,8 @@ namespace Step32
 
     std::vector<double> rhs_values(n_q_points);
 
-    LA::MPI::Vector
-    rhs (temperature_mass_matrix.row_partitioner()),
-        solution (temperature_mass_matrix.row_partitioner());
+    LA::MPI::Vector rhs (temperature_dof_handler.locally_owned_dofs(), MPI_COMM_WORLD);
+    LA::MPI::Vector solution (temperature_dof_handler.locally_owned_dofs(), MPI_COMM_WORLD);
 
     const EquationData::TemperatureInitialValues<dim> initial_temperature;
 
@@ -1763,7 +1762,11 @@ namespace Step32
     SolverCG<LA::MPI::Vector> cg(solver_control);
 
     LA::MPI::PreconditionJacobi preconditioner_mass;
+#ifdef USE_PETSC_LA
+    preconditioner_mass.initialize(temperature_mass_matrix);
+#else
     preconditioner_mass.initialize(temperature_mass_matrix, 1.3);
+#endif
 
     cg.solve (temperature_mass_matrix, solution, rhs, preconditioner_mass);
 

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.