]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Cleanup AdditionalData in SolverCG and SolverGMRES 4510/head
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sat, 10 Jun 2017 09:48:35 +0000 (11:48 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Tue, 13 Jun 2017 09:17:29 +0000 (11:17 +0200)
doc/news/changes/incompatibilities/20170612DanielArndt-1 [new file with mode: 0644]
include/deal.II/lac/solver_cg.h
include/deal.II/lac/solver_gmres.h
tests/lac/gmres_eigenvalues.cc
tests/lac/gmres_eigenvalues.output
tests/lac/solver_signals.output
tests/lapack/solver_cg.cc
tests/lapack/solver_cg.output
tests/matrix_free/estimate_condition_number_mass.cc

diff --git a/doc/news/changes/incompatibilities/20170612DanielArndt-1 b/doc/news/changes/incompatibilities/20170612DanielArndt-1
new file mode 100644 (file)
index 0000000..f8e8b64
--- /dev/null
@@ -0,0 +1,5 @@
+Changed: The deprecated data in SolverCG::AdditionalData
+and SolverGMRES::AdditionalData::compute_eigenvalues have been removed.
+Use the respective connect_* member functions instead.
+<br>
+(Daniel Arndt, 2017/06/12)
index 53fcf63bdfabe84ac43149b14de70b63225a9772..87a690154bcfd1373d96443625b5a7a50e077f82 100644 (file)
@@ -45,11 +45,6 @@ class PreconditionIdentity;
  * solution vector must be passed as template argument, and defaults to
  * dealii::Vector<double>.
  *
- * Like all other solver classes, this class has a local structure called @p
- * AdditionalData which is used to pass additional parameters to the solver.
- * For this class, there is (among other things) a switch allowing for
- * additional output for the computation of eigenvalues of the matrix.
- *
  * @note This version of CG is taken from D. Braess's book "Finite Elements".
  * It requires a symmetric preconditioner (i.e., for example, SOR is not a
  * possible choice).
@@ -83,9 +78,6 @@ class PreconditionIdentity;
  * connect_condition_number_slot. These slots will then be called from the
  * solver with the estimates as argument.
  *
- * @deprecated Alternatively these estimates can be written to deallog by
- * setting flags in @p AdditionalData.
- *
  * <h3>Observing the progress of linear solver iterations</h3>
  *
  * The solve() function of this class uses the mechanism described in the
@@ -106,52 +98,10 @@ public:
 
   /**
    * Standardized data struct to pipe additional data to the solver.
+   * Here, it doesn't store anything but just exists for consistency
+   * with the other solver classes.
    */
-  struct AdditionalData
-  {
-    /**
-     * Write coefficients alpha and beta to the log file for later use in
-     * eigenvalue estimates.
-     */
-    bool log_coefficients;
-
-    /**
-     * Compute the condition number of the projected matrix.
-     *
-     * @note Requires LAPACK support.
-     */
-    bool compute_condition_number;
-
-    /**
-     * Compute the condition number of the projected matrix in each step.
-     *
-     * @note Requires LAPACK support.
-     */
-    bool compute_all_condition_numbers;
-
-    /**
-     * Compute all eigenvalues of the projected matrix.
-     *
-     * @note Requires LAPACK support.
-     */
-    bool compute_eigenvalues;
-
-    /**
-     * Constructor. Initialize data fields.  Confer the description of those.
-     * @deprecated Instead use: connect_coefficients_slot,
-     * connect_condition_number_slot, and connect_eigenvalues_slot.
-     */
-    explicit
-    AdditionalData (const bool log_coefficients,
-                    const bool compute_condition_number = false,
-                    const bool compute_all_condition_numbers = false,
-                    const bool compute_eigenvalues = false) DEAL_II_DEPRECATED;
-
-    /**
-     * Constructor. Initializes all data fields to false.
-     */
-    AdditionalData();
-  };
+  struct AdditionalData {};
 
   /**
    * Constructor.
@@ -228,17 +178,13 @@ protected:
    * Estimates the eigenvalues from diagonal and offdiagonal. Uses these
    * estimate to compute the condition number. Calls the signals
    * eigenvalues_signal and cond_signal with these estimates as arguments.
-   * Outputs the eigenvalues/condition-number to deallog if
-   * log_eigenvalues/log_cond is true.
    */
   static void
   compute_eigs_and_cond(
     const std::vector<double> &diagonal,
     const std::vector<double> &offdiagonal,
     const boost::signals2::signal<void (const std::vector<double> &)> &eigenvalues_signal,
-    const boost::signals2::signal<void (double)> &cond_signal,
-    const bool log_eigenvalues,
-    const bool log_cond);
+    const boost::signals2::signal<void (double)> &cond_signal);
 
   /**
    * Temporary vectors, allocated through the @p VectorMemory object at the
@@ -292,35 +238,6 @@ private:
 
 #ifndef DOXYGEN
 
-template <typename VectorType>
-inline
-SolverCG<VectorType>::AdditionalData::
-AdditionalData (const bool log_coefficients,
-                const bool compute_condition_number,
-                const bool compute_all_condition_numbers,
-                const bool compute_eigenvalues)
-  :
-  log_coefficients (log_coefficients),
-  compute_condition_number(compute_condition_number),
-  compute_all_condition_numbers(compute_all_condition_numbers),
-  compute_eigenvalues(compute_eigenvalues)
-{}
-
-
-
-template <typename VectorType>
-inline
-SolverCG<VectorType>::AdditionalData::
-AdditionalData ()
-  :
-  log_coefficients (false),
-  compute_condition_number(false),
-  compute_all_condition_numbers(false),
-  compute_eigenvalues(false)
-{}
-
-
-
 template <typename VectorType>
 SolverCG<VectorType>::SolverCG (SolverControl        &cn,
                                 VectorMemory<VectorType> &mem,
@@ -382,13 +299,10 @@ SolverCG<VectorType>::compute_eigs_and_cond
 (const std::vector<double> &diagonal,
  const std::vector<double> &offdiagonal,
  const boost::signals2::signal<void (const std::vector<double> &)> &eigenvalues_signal,
- const boost::signals2::signal<void (double)>                      &cond_signal,
- const bool                log_eigenvalues,
- const bool                log_cond)
+ const boost::signals2::signal<void (double)>                      &cond_signal)
 {
   //Avoid computing eigenvalues unless they are needed.
-  if (!cond_signal.empty()|| !eigenvalues_signal.empty()  || log_cond ||
-      log_eigenvalues)
+  if (!cond_signal.empty()|| !eigenvalues_signal.empty())
     {
       TridiagonalMatrix<double> T(diagonal.size(), true);
       for (size_type i=0; i<diagonal.size(); ++i)
@@ -403,12 +317,6 @@ SolverCG<VectorType>::compute_eigs_and_cond
         {
           double condition_number=T.eigenvalue(T.n()-1)/T.eigenvalue(0);
           cond_signal(condition_number);
-          //Send to deallog
-          if (log_cond)
-            {
-              deallog << "Condition number estimate: " <<
-                      condition_number << std::endl;
-            }
         }
       //Avoid copying the eigenvalues of T to a vector unless a signal is
       //connected.
@@ -421,12 +329,6 @@ SolverCG<VectorType>::compute_eigs_and_cond
             }
           eigenvalues_signal(eigenvalues);
         }
-      if (log_eigenvalues)
-        {
-          for (size_type i=0; i<T.n(); ++i)
-            deallog << ' ' << T.eigenvalue(i);
-          deallog << std::endl;
-        }
     }
 
 }
@@ -454,10 +356,7 @@ SolverCG<VectorType>::solve (const MatrixType         &A,
   const bool do_eigenvalues = !condition_number_signal.empty()
                               ||!all_condition_numbers_signal.empty()
                               ||!eigenvalues_signal.empty()
-                              ||!all_eigenvalues_signal.empty()
-                              || additional_data.compute_condition_number
-                              || additional_data.compute_all_condition_numbers
-                              || additional_data.compute_eigenvalues;
+                              ||!all_eigenvalues_signal.empty();
 
   // vectors used for eigenvalue
   // computations
@@ -555,8 +454,6 @@ SolverCG<VectorType>::solve (const MatrixType         &A,
             }
 
           this->coefficients_signal(alpha,beta);
-          if (additional_data.log_coefficients)
-            deallog << "alpha-beta:" << alpha << '\t' << beta << std::endl;
           // set up the vectors
           // containing the diagonal
           // and the off diagonal of
@@ -568,8 +465,7 @@ SolverCG<VectorType>::solve (const MatrixType         &A,
               offdiagonal.push_back(std::sqrt(beta)/alpha);
             }
           compute_eigs_and_cond(diagonal,offdiagonal,all_eigenvalues_signal,
-                                all_condition_numbers_signal,false,
-                                additional_data.compute_all_condition_numbers);
+                                all_condition_numbers_signal);
         }
     }
   catch (...)
@@ -578,10 +474,7 @@ SolverCG<VectorType>::solve (const MatrixType         &A,
       throw;
     }
   compute_eigs_and_cond(diagonal,offdiagonal,eigenvalues_signal,
-                        condition_number_signal,
-                        additional_data.compute_eigenvalues,
-                        (additional_data.compute_condition_number &&
-                         !additional_data.compute_all_condition_numbers));
+                        condition_number_signal);
 
   // Deallocate Memory
   cleanup();
index bd15eb78ee93b071c7049be67471639d89e5d34d..87b9118bce840b285c528c5e033514d08617aa58 100644 (file)
@@ -198,17 +198,6 @@ public:
                     const bool use_default_residual = true,
                     const bool force_re_orthogonalization = false);
 
-    /**
-     * Constructor.
-     * @deprecated To obtain the estimated eigenvalues instead use:
-     * connect_eigenvalues_slot
-     */
-    AdditionalData (const unsigned int max_n_tmp_vectors,
-                    const bool right_preconditioning,
-                    const bool use_default_residual,
-                    const bool force_re_orthogonalization,
-                    const bool compute_eigenvalues) DEAL_II_DEPRECATED;
-
     /**
      * Maximum number of temporary vectors. This parameter controls the size
      * of the Arnoldi basis, which for historical reasons is
@@ -237,17 +226,6 @@ public:
      * if necessary.
      */
     bool force_re_orthogonalization;
-
-    /**
-     * Compute all eigenvalues of the Hessenberg matrix generated while
-     * solving, i.e., the projected system matrix. This gives an approximation
-     * of the eigenvalues of the (preconditioned) system matrix. Since the
-     * Hessenberg matrix is thrown away at restart, the eigenvalues are
-     * printed for every 30 iterations.
-     *
-     * @note Requires LAPACK support.
-     */
-    bool compute_eigenvalues;
   };
 
   /**
@@ -407,8 +385,7 @@ protected:
    * Estimates the eigenvalues from the Hessenberg matrix, H_orig, generated
    * during the inner iterations. Uses these estimate to compute the condition
    * number. Calls the signals eigenvalues_signal and cond_signal with these
-   * estimates as arguments. Outputs the eigenvalues to deallog if
-   * log_eigenvalues is true.
+   * estimates as arguments.
    */
   static void
   compute_eigs_and_cond(
@@ -416,8 +393,7 @@ protected:
     const unsigned int dim,
     const boost::signals2::signal<void (const std::vector<std::complex<double> > &)> &eigenvalues_signal,
     const boost::signals2::signal<void (const FullMatrix<double> &)> &hessenberg_signal,
-    const boost::signals2::signal<void(double)> &cond_signal,
-    const bool log_eigenvalues);
+    const boost::signals2::signal<void(double)> &cond_signal);
 
   /**
    * Projected system matrix
@@ -616,26 +592,7 @@ AdditionalData (const unsigned int max_n_tmp_vectors,
   max_n_tmp_vectors(max_n_tmp_vectors),
   right_preconditioning(right_preconditioning),
   use_default_residual(use_default_residual),
-  force_re_orthogonalization(force_re_orthogonalization),
-  compute_eigenvalues(false)
-{}
-
-
-
-template <class VectorType>
-inline
-SolverGMRES<VectorType>::AdditionalData::
-AdditionalData (const unsigned int max_n_tmp_vectors,
-                const bool         right_preconditioning,
-                const bool         use_default_residual,
-                const bool         force_re_orthogonalization,
-                const bool         compute_eigenvalues)
-  :
-  max_n_tmp_vectors(max_n_tmp_vectors),
-  right_preconditioning(right_preconditioning),
-  use_default_residual(use_default_residual),
-  force_re_orthogonalization(force_re_orthogonalization),
-  compute_eigenvalues(compute_eigenvalues)
+  force_re_orthogonalization(force_re_orthogonalization)
 {}
 
 
@@ -758,12 +715,11 @@ SolverGMRES<VectorType>::compute_eigs_and_cond
  const unsigned int                           dim,
  const boost::signals2::signal<void (const std::vector<std::complex<double> > &)> &eigenvalues_signal,
  const boost::signals2::signal<void (const FullMatrix<double> &)> &hessenberg_signal,
- const boost::signals2::signal<void (double)> &cond_signal,
- const bool                                   log_eigenvalues)
+ const boost::signals2::signal<void (double)> &cond_signal)
 {
   //Avoid copying the Hessenberg matrix if it isn't needed.
   if (!eigenvalues_signal.empty() || !hessenberg_signal.empty()
-      || !cond_signal.empty() || log_eigenvalues )
+      || !cond_signal.empty())
     {
       LAPACKFullMatrix<double> mat(dim,dim);
       for (unsigned int i=0; i<dim; ++i)
@@ -771,7 +727,7 @@ SolverGMRES<VectorType>::compute_eigs_and_cond
           mat(i,j) = H_orig(i,j);
       hessenberg_signal(H_orig);
       //Avoid computing eigenvalues if they are not needed.
-      if (!eigenvalues_signal.empty() || log_eigenvalues )
+      if (!eigenvalues_signal.empty())
         {
           //Copy mat so that we can compute svd below. Necessary since
           //compute_eigenvalues will leave mat in state LAPACKSupport::unusable.
@@ -784,13 +740,6 @@ SolverGMRES<VectorType>::compute_eigs_and_cond
           std::sort(eigenvalues.begin(), eigenvalues.end(),
                     internal::SolverGMRES::complex_less_pred);
           eigenvalues_signal(eigenvalues);
-          if (log_eigenvalues)
-            {
-              deallog << "Eigenvalue estimate: ";
-              for (unsigned int i=0; i<mat_eig.n(); ++i)
-                deallog << ' ' << eigenvalues[i];
-              deallog << std::endl;
-            }
         }
       //Calculate condition number, avoid calculating the svd if a slot
       //isn't connected. Need at least a 2-by-2 matrix to do the estimate.
@@ -837,8 +786,7 @@ SolverGMRES<VectorType>::solve (const MatrixType         &A,
     ||!eigenvalues_signal.empty()
     ||!all_eigenvalues_signal.empty()
     ||!hessenberg_signal.empty()
-    ||!all_hessenberg_signal.empty()
-    ||additional_data.compute_eigenvalues;
+    ||!all_hessenberg_signal.empty();
   // for eigenvalue computation, need to collect the Hessenberg matrix (before
   // applying Givens rotations)
   FullMatrix<double> H_orig;
@@ -1064,9 +1012,7 @@ SolverGMRES<VectorType>::solve (const MatrixType         &A,
           H1(i,j) = H(i,j);
 
       compute_eigs_and_cond(H_orig,dim,all_eigenvalues_signal,
-                            all_hessenberg_signal,
-                            all_condition_numbers_signal,
-                            additional_data.compute_eigenvalues);
+                            all_hessenberg_signal, condition_number_signal);
 
       H1.backward(h,gamma);
 
@@ -1087,8 +1033,7 @@ SolverGMRES<VectorType>::solve (const MatrixType         &A,
   while (iteration_state == SolverControl::iterate);
 
   compute_eigs_and_cond(H_orig,dim,eigenvalues_signal,hessenberg_signal,
-                        condition_number_signal,
-                        false);
+                        condition_number_signal);
 
   if (!krylov_space_signal.empty())
     krylov_space_signal(tmp_vectors);
index 9c12bb670628e61bca66dc08466017fc7892585e..6c29b79706fe61449c7dae168bc9fa5a49f28c32 100644 (file)
 #include <deal.II/lac/solver_gmres.h>
 #include <deal.II/lac/precondition.h>
 
-
+template<class NUMBER>
+void output_eigenvalues(const std::vector<NUMBER> &eigenvalues,const std::string &text)
+{
+  deallog<< text;
+  for (unsigned int j = 0; j < eigenvalues.size(); ++j)
+    {
+      deallog<< ' ' << eigenvalues.at(j);
+    }
+  deallog << std::endl;
+}
 
 template <typename number>
 void test (unsigned int variant)
@@ -66,16 +75,19 @@ void test (unsigned int variant)
   SolverControl control(1000, variant==1?1e-4:1e-13);
   typename SolverGMRES<Vector<number> >::AdditionalData data;
   data.max_n_tmp_vectors = 80;
-  data.compute_eigenvalues = true;
 
   SolverGMRES<Vector<number> > solver(control, data);
+  solver.connect_eigenvalues_slot(
+    std::bind(output_eigenvalues<std::complex<double>>,
+              std::placeholders::_1,"Eigenvalue estimate: "));
   solver.solve(matrix, sol, rhs, PreconditionIdentity());
 
   if (variant == 0)
     {
-      typename SolverCG<Vector<number> >::AdditionalData cg_data;
-      cg_data.compute_eigenvalues = true;
-      SolverCG<Vector<number> > solver_cg(control, cg_data);
+      SolverCG<Vector<number> > solver_cg(control);
+      solver_cg.connect_eigenvalues_slot(
+        std::bind(output_eigenvalues<double>,
+                  std::placeholders::_1,"Eigenvalue estimate: "));
       sol = 0;
       solver_cg.solve(matrix, sol, rhs, PreconditionIdentity());
     }
index 4bcadb4367346c27e03036298cbb53b456a19b4f..1b6aa8eb44f2bbc3a690fef48c3e12ca34232fa6 100644 (file)
@@ -4,7 +4,7 @@ DEAL:double:0:GMRES::Convergence step 56 value 0
 DEAL:double:0:GMRES::Eigenvalue estimate:  (1.0,0.0) (2.0,0.0) (3.0,0.0) (4.0,0.0) (5.0,0.0) (6.0,0.0) (7.0,0.0) (8.0,0.0) (9.0,0.0) (10.,0.0) (11.,0.0) (12.,0.0) (13.,0.0) (14.,0.0) (15.,0.0) (16.,0.0) (17.,0.0) (18.,0.0) (19.,0.0) (21.,0.0) (22.,0.0) (23.,0.0) (25.,0.0) (26.,0.0) (27.,0.0) (29.,0.0) (30.,0.0) (32.,0.0) (33.,0.0) (35.,0.0) (36.,0.0) (38.,0.0) (39.,0.0) (40.,0.0) (42.,0.0) (43.,0.0) (44.,0.0) (46.,0.0) (47.,0.0) (48.,0.0) (49.,0.0) (50.,0.0) (51.,0.0) (52.,0.0) (53.,0.0) (54.,0.0) (55.,0.0) (56.,0.0) (57.,0.0) (58.,0.0) (59.,0.0) (60.,0.0) (61.,0.0) (62.,0.0) (63.,0.0) (64.,0.0)
 DEAL:double:0:cg::Starting value 8.0
 DEAL:double:0:cg::Convergence step 56 value 0
-DEAL:double:0:cg:: 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10. 11. 12. 13. 14. 15. 16. 17. 18. 20. 21. 22. 24. 25. 27. 28. 30. 31. 33. 34. 35. 37. 38. 40. 41. 43. 44. 45. 47. 48. 49. 50. 51. 52. 53. 54. 55. 56. 57. 58. 59. 60. 61. 62. 63. 64.
+DEAL:double:0:cg::Eigenvalue estimate:  1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10. 11. 12. 13. 14. 15. 16. 17. 18. 20. 21. 22. 24. 25. 27. 28. 30. 31. 33. 34. 35. 37. 38. 40. 41. 43. 44. 45. 47. 48. 49. 50. 51. 52. 53. 54. 55. 56. 57. 58. 59. 60. 61. 62. 63. 64.
 DEAL:double:1:GMRES::Starting value 4.0
 DEAL:double:1:GMRES::Convergence step 16 value 0
 DEAL:double:1:GMRES::Eigenvalue estimate:  (1.0,0.0) (16.,0.0) (81.,0.0) (2.6e+02,0.0) (6.3e+02,0.0) (1.3e+03,0.0) (2.4e+03,0.0) (4.1e+03,0.0) (6.6e+03,0.0) (1.0e+04,0.0) (1.5e+04,0.0) (2.1e+04,0.0) (2.9e+04,0.0) (3.8e+04,0.0) (5.1e+04,0.0) (6.6e+04,0.0)
index 896dae47323ebde8261f6ee24ab9e94a27e6cf7a..f597ef56d0922d0156cf9e64bfd0c16e147167f3 100644 (file)
@@ -121,12 +121,12 @@ DEAL:cg::Final condition number: 364.1
 DEAL:cg::Final Eigenvalues:  0.02191 0.1088 0.1958 0.2789 0.3660 0.5253 0.6175 0.8323 0.9462 1.192 1.322 1.598 1.804 2.024 2.402 2.518 2.834 3.186 3.485 3.728 3.962 4.346 4.672 4.910 5.218 5.562 5.863 6.072 6.333 6.629 6.870 7.024 7.222 7.427 7.590 7.663 7.801 7.891 7.978
 DEAL:GMRES::Starting value 29.00
 DEAL:GMRES::Eigenvalues:  (0.02191,0.000) (0.1089,0.000) (0.2179,0.000) (0.2902,0.000) (0.5179,0.000) (0.6680,0.000) (0.8533,0.000) (1.202,0.000) (1.587,0.000) (1.962,0.000) (2.224,0.000) (2.520,0.000) (2.908,0.000) (3.304,0.000) (3.685,0.000) (4.034,0.000) (4.554,0.000) (5.012,0.000) (5.447,0.000) (5.860,0.000) (6.248,0.000) (6.609,0.000) (6.934,0.000) (7.216,0.000) (7.456,0.000) (7.665,0.000) (7.836,0.000) (7.954,0.000)
-DEAL:GMRES::Condition number: 363.0
+DEAL:GMRES::Final condition number: 363.0
 DEAL:GMRES::Eigenvalues:  (0.02765,0.000) (0.1418,0.000) (0.2014,0.000) (0.3484,0.000) (0.5334,0.000) (0.6447,0.000) (1.032,0.000) (1.319,0.000) (1.634,0.000) (1.937,0.000) (2.310,0.000) (2.760,0.000) (3.205,0.000) (3.694,0.000) (4.075,0.000) (4.475,0.000) (4.910,0.000) (5.333,0.000) (5.741,0.000) (6.146,0.000) (6.553,0.000) (6.901,0.000) (7.129,0.000) (7.393,0.000) (7.647,0.000) (7.797,0.000) (7.891,0.000) (7.978,0.000)
-DEAL:GMRES::Condition number: 288.5
+DEAL:GMRES::Final condition number: 288.5
 DEAL:GMRES::Convergence step 65 value 0.0009799
 DEAL:GMRES::Eigenvalues:  (0.02681,0.000) (0.2333,0.000) (1.178,0.000) (2.161,0.000) (3.554,0.000) (4.930,0.000) (6.307,0.000) (7.170,0.000) (7.726,0.000)
-DEAL:GMRES::Condition number: 288.2
+DEAL:GMRES::Final condition number: 288.2
 DEAL:GMRES::Hessenberg matrix: 
 DEAL:GMRES::0.1982 0.8628 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
 DEAL:GMRES::0.8628 6.098 1.971 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
index f4cc0337c252147839cc8a524e9048969830c054..d84dc6e1195a2ca1697e1894286c74a5e25b4b66 100644 (file)
 #include <deal.II/lac/solver_cg.h>
 #include <deal.II/lac/precondition.h>
 
+void output_double_number(double input,const std::string &text)
+{
+  deallog<<text<< input<<std::endl;
+}
+
+template<class NUMBER>
+void output_eigenvalues(const std::vector<NUMBER> &eigenvalues,const std::string &text)
+{
+  deallog<< text;
+  for (unsigned int j = 0; j < eigenvalues.size(); ++j)
+    {
+      deallog<< ' ' << eigenvalues.at(j);
+    }
+  deallog << std::endl;
+}
+
+
 template<typename SolverType, typename MatrixType, typename VectorType, class PRECONDITION>
 void
 check_solve (SolverType         &solver,
@@ -82,8 +99,12 @@ int main()
   GrowingVectorMemory<> mem;
   SolverControl control(100, 1.e-3);
   SolverControl verbose_control(100, 1.e-3, true);
-  SolverCG<>::AdditionalData data(false, true, true, true);
-  SolverCG<> cg(control, mem, data);
+  SolverCG<> cg(control, mem);
+  cg.connect_condition_number_slot(
+    std::bind(output_double_number,std::placeholders::_1,"Condition number estimate: "),true);
+  cg.connect_eigenvalues_slot(
+    std::bind(output_eigenvalues<double>,std::placeholders::_1,"Final Eigenvalues: "));
+
 
   for (unsigned int size=4; size <= 30; size *= 3)
     {
index 609ea6e9a6f28a6ba3e04de49059e7ef1f590eb9..69e9273110c3795695d0a7f5b504daa5de5d3686 100644 (file)
@@ -3,20 +3,20 @@ DEAL::Size 4 Unknowns 9
 DEAL:no-fail:cg::Starting value 3.000
 DEAL:no-fail:cg::Condition number estimate: 3.535
 DEAL:no-fail:cg::Convergence step 3 value 0
-DEAL:no-fail:cg:: 1.176 4.157
+DEAL:no-fail:cg::Final Eigenvalues:  1.176 4.157
 DEAL:no:cg::Starting value 3.000
 DEAL:no:cg::Condition number estimate: 3.535
 DEAL:no:cg::Convergence step 3 value 0
-DEAL:no:cg:: 1.176 4.157
+DEAL:no:cg::Final Eigenvalues:  1.176 4.157
 DEAL:rich:cg::Starting value 3.000
 DEAL:rich:cg::Condition number estimate: 3.535
 DEAL:rich:cg::Convergence step 3 value 0
-DEAL:rich:cg:: 0.7056 2.494
+DEAL:rich:cg::Final Eigenvalues:  0.7056 2.494
 DEAL:ssor:cg::Starting value 3.000
 DEAL:ssor:cg::Condition number estimate: 1.384
 DEAL:ssor:cg::Condition number estimate: 1.421
 DEAL:ssor:cg::Convergence step 4 value 6.018e-05
-DEAL:ssor:cg:: 0.7030 0.8172 0.9988
+DEAL:ssor:cg::Final Eigenvalues:  0.7030 0.8172 0.9988
 DEAL::Size 12 Unknowns 121
 DEAL:no-fail:cg::Starting value 11.00
 DEAL:no-fail:cg::Condition number estimate: 11.01
@@ -28,7 +28,7 @@ DEAL:no-fail:cg::Condition number estimate: 51.74
 DEAL:no-fail:cg::Condition number estimate: 53.54
 DEAL:no-fail:cg::Condition number estimate: 54.75
 DEAL:no-fail:cg::Failure step 10 value 0.1496
-DEAL:no-fail:cg:: 0.1363 0.6565 1.499 2.357 3.131 3.994 5.392 6.576 7.464
+DEAL:no-fail:cg::Final Eigenvalues:  0.1363 0.6565 1.499 2.357 3.131 3.994 5.392 6.576 7.464
 DEAL:no-fail::Exception: SolverControl::NoConvergence (it, res)
 DEAL:no:cg::Starting value 11.00
 DEAL:no:cg::Condition number estimate: 11.01
@@ -45,7 +45,7 @@ DEAL:no:cg::Condition number estimate: 57.32
 DEAL:no:cg::Condition number estimate: 57.64
 DEAL:no:cg::Condition number estimate: 57.68
 DEAL:no:cg::Convergence step 15 value 0.0005794
-DEAL:no:cg:: 0.1363 0.6539 1.173 1.552 2.119 2.616 3.388 3.933 4.821 5.350 6.003 6.785 7.337 7.862
+DEAL:no:cg::Final Eigenvalues:  0.1363 0.6539 1.173 1.552 2.119 2.616 3.388 3.933 4.821 5.350 6.003 6.785 7.337 7.862
 DEAL:rich:cg::Starting value 11.00
 DEAL:rich:cg::Condition number estimate: 11.01
 DEAL:rich:cg::Condition number estimate: 21.10
@@ -61,7 +61,7 @@ DEAL:rich:cg::Condition number estimate: 57.32
 DEAL:rich:cg::Condition number estimate: 57.64
 DEAL:rich:cg::Condition number estimate: 57.68
 DEAL:rich:cg::Convergence step 15 value 0.0005794
-DEAL:rich:cg:: 0.08178 0.3924 0.7039 0.9313 1.271 1.569 2.033 2.360 2.893 3.210 3.602 4.071 4.402 4.717
+DEAL:rich:cg::Final Eigenvalues:  0.08178 0.3924 0.7039 0.9313 1.271 1.569 2.033 2.360 2.893 3.210 3.602 4.071 4.402 4.717
 DEAL:ssor:cg::Starting value 11.00
 DEAL:ssor:cg::Condition number estimate: 4.354
 DEAL:ssor:cg::Condition number estimate: 5.511
@@ -70,4 +70,4 @@ DEAL:ssor:cg::Condition number estimate: 5.706
 DEAL:ssor:cg::Condition number estimate: 5.715
 DEAL:ssor:cg::Condition number estimate: 5.718
 DEAL:ssor:cg::Convergence step 8 value 0.0005816
-DEAL:ssor:cg:: 0.1748 0.4549 0.5519 0.6780 0.7983 0.9058 0.9996
+DEAL:ssor:cg::Final Eigenvalues:  0.1748 0.4549 0.5519 0.6780 0.7983 0.9058 0.9996
index d8e53ab82203a6417301bf96e7927eddda276717..8cc8f6472a7d6903f47108b92e14f4698346fed4 100644 (file)
@@ -37,6 +37,11 @@ std::ofstream logfile("output");
 #include <deal.II/base/function_lib.h>
 
 
+void output_double_number(double input,const std::string &text)
+{
+  deallog<<text<< input<<std::endl;
+}
+
 template <int dim, int fe_degree, typename Number>
 void
 mass_operator (const MatrixFree<dim,Number>  &data,
@@ -127,9 +132,9 @@ void test (const FiniteElement<dim> &fe,
   // accumulate differently. Beware of this strange solver setting when seeing
   // "failure" in the output
   SolverControl control(n_iterations, 0);
-  typename SolverCG<>::AdditionalData data;
-  data.compute_condition_number = true;
-  SolverCG<> solver(control,data);
+  SolverCG<> solver(control);
+  solver.connect_condition_number_slot(
+    std::bind(output_double_number,std::placeholders::_1,"Condition number estimate: "));
   try
     {
       solver.solve(mf, out, in, PreconditionIdentity());

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.