]> https://gitweb.dealii.org/ - dealii.git/commitdiff
extend mpi/step-40 to complex-valued PETSc
authorDenis Davydov <davydden@gmail.com>
Mon, 28 Mar 2016 20:41:55 +0000 (22:41 +0200)
committerDenis Davydov <davydden@gmail.com>
Tue, 29 Mar 2016 18:48:46 +0000 (20:48 +0200)
tests/mpi/step-40.cc
tests/mpi/step-40.with_petsc=true.mpirun=10.output.2 [new file with mode: 0644]
tests/mpi/step-40.with_petsc=true.mpirun=3.output.2 [new file with mode: 0644]
tests/mpi/step-40.with_petsc=true.mpirun=4.output.2 [new file with mode: 0644]

index b8927a4f76aac2b01f57748e2a484f538fa90444..619df472be458bac1597839e0d20521aa3ac4ab8 100644 (file)
@@ -138,7 +138,7 @@ namespace Step40
                                       locally_relevant_dofs,
                                       mpi_communicator);
     system_rhs.reinit (locally_owned_dofs, mpi_communicator);
-    system_rhs = 0;
+    system_rhs = PetscScalar();
 
     constraints.clear ();
     constraints.reinit (locally_relevant_dofs);
@@ -182,8 +182,8 @@ namespace Step40
     const unsigned int   dofs_per_cell = fe.dofs_per_cell;
     const unsigned int   n_q_points    = quadrature_formula.size();
 
-    FullMatrix<double>   cell_matrix (dofs_per_cell, dofs_per_cell);
-    Vector<double>       cell_rhs (dofs_per_cell);
+    FullMatrix<PetscScalar>   cell_matrix (dofs_per_cell, dofs_per_cell);
+    Vector<PetscScalar>       cell_rhs (dofs_per_cell);
 
     std::vector<types::global_dof_index> local_dof_indices (dofs_per_cell);
 
@@ -193,8 +193,8 @@ namespace Step40
     for (; cell!=endc; ++cell)
       if (cell->is_locally_owned())
         {
-          cell_matrix = 0;
-          cell_rhs = 0;
+          cell_matrix = PetscScalar();
+          cell_rhs = PetscScalar();
 
           fe_values.reinit (cell);
 
@@ -248,6 +248,7 @@ namespace Step40
 
     PETScWrappers::SolverCG solver(solver_control, mpi_communicator);
 
+#ifndef PETSC_USE_COMPLEX
     PETScWrappers::PreconditionBoomerAMG
     preconditioner(system_matrix,
                    PETScWrappers::PreconditionBoomerAMG::AdditionalData(true));
@@ -256,6 +257,12 @@ namespace Step40
       solver.solve (system_matrix, completely_distributed_solution, system_rhs,
                     preconditioner),
       solver_control.last_step(), 10, 10);
+#else
+    check_solver_within_range(
+      solver.solve (system_matrix, completely_distributed_solution, system_rhs,
+                    PETScWrappers::PreconditionJacobi(system_matrix)),
+      solver_control.last_step(), 120, 260);
+#endif
 
     pcout << "   Solved in " << solver_control.last_step()
           << " iterations." << std::endl;
diff --git a/tests/mpi/step-40.with_petsc=true.mpirun=10.output.2 b/tests/mpi/step-40.with_petsc=true.mpirun=10.output.2
new file mode 100644 (file)
index 0000000..c71dce0
--- /dev/null
@@ -0,0 +1,17 @@
+
+Cycle 0:
+   Number of active cells:       1024
+      104+100+104+100+104+104+100+104+100+104+
+   Number of degrees of freedom: 4225
+      465+416+432+416+416+432+416+416+400+416+
+DEAL:mpi::Solver stopped within 120 - 260 iterations
+   Solved in 126 iterations.
+
+Cycle 1:
+   Number of active cells:       4096
+      408+412+408+412+408+408+412+408+412+408+
+   Number of degrees of freedom: 16641
+      1729+1680+1664+1680+1632+1664+1680+1632+1648+1632+
+DEAL:mpi::Solver stopped within 120 - 260 iterations
+   Solved in 251 iterations.
+
diff --git a/tests/mpi/step-40.with_petsc=true.mpirun=3.output.2 b/tests/mpi/step-40.with_petsc=true.mpirun=3.output.2
new file mode 100644 (file)
index 0000000..ea0c6a3
--- /dev/null
@@ -0,0 +1,17 @@
+
+Cycle 0:
+   Number of active cells:       1024
+      340+344+340+
+   Number of degrees of freedom: 4225
+      1453+1412+1360+
+DEAL:mpi::Solver stopped within 120 - 260 iterations
+   Solved in 126 iterations.
+
+Cycle 1:
+   Number of active cells:       4096
+      1364+1368+1364+
+   Number of degrees of freedom: 16641
+      5645+5540+5456+
+DEAL:mpi::Solver stopped within 120 - 260 iterations
+   Solved in 251 iterations.
+
diff --git a/tests/mpi/step-40.with_petsc=true.mpirun=4.output.2 b/tests/mpi/step-40.with_petsc=true.mpirun=4.output.2
new file mode 100644 (file)
index 0000000..41f1432
--- /dev/null
@@ -0,0 +1,17 @@
+
+Cycle 0:
+   Number of active cells:       1024
+      256+256+256+256+
+   Number of degrees of freedom: 4225
+      1089+1056+1056+1024+
+DEAL:mpi::Solver stopped within 120 - 260 iterations
+   Solved in 126 iterations.
+
+Cycle 1:
+   Number of active cells:       4096
+      1024+1024+1024+1024+
+   Number of degrees of freedom: 16641
+      4225+4160+4160+4096+
+DEAL:mpi::Solver stopped within 120 - 260 iterations
+   Solved in 251 iterations.
+

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.