From: Wolfgang Bangerth Date: Thu, 5 Oct 2017 22:57:41 +0000 (-0600) Subject: Replace logstream.push/pop by Logstream::Prefix. X-Git-Tag: v9.0.0-rc1~973^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4931cc6432e52a6ab81911d386b0e66cc09c5465;p=dealii.git Replace logstream.push/pop by Logstream::Prefix. --- diff --git a/include/deal.II/lac/solver_bicgstab.h b/include/deal.II/lac/solver_bicgstab.h index 7450a87142..43093e076d 100644 --- a/include/deal.II/lac/solver_bicgstab.h +++ b/include/deal.II/lac/solver_bicgstab.h @@ -449,7 +449,7 @@ SolverBicgstab::solve(const MatrixType &A, const VectorType &b, const PreconditionerType &preconditioner) { - deallog.push("Bicgstab"); + LogStream::Prefix prefix("Bicgstab"); Vr = typename VectorMemory::Pointer(this->memory); Vrbar = typename VectorMemory::Pointer(this->memory); Vp = typename VectorMemory::Pointer(this->memory); @@ -488,8 +488,6 @@ SolverBicgstab::solve(const MatrixType &A, } while (state.breakdown == true); - deallog.pop(); - // in case of failure: throw exception AssertThrow(state.state == SolverControl::success, SolverControl::NoConvergence (state.last_step, diff --git a/include/deal.II/lac/solver_cg.h b/include/deal.II/lac/solver_cg.h index 3d784e5867..028d4031e5 100644 --- a/include/deal.II/lac/solver_cg.h +++ b/include/deal.II/lac/solver_cg.h @@ -310,7 +310,7 @@ SolverCG::solve (const MatrixType &A, { SolverControl::State conv=SolverControl::iterate; - deallog.push("cg"); + LogStream::Prefix prefix("cg"); // Memory allocation typename VectorMemory::Pointer g_pointer(this->memory); @@ -336,114 +336,102 @@ SolverCG::solve (const MatrixType &A, int it=0; double res = -std::numeric_limits::max(); - try + double eigen_beta_alpha = 0; + + // resize the vectors, but do not set + // the values since they'd be overwritten + // soon anyway. + g.reinit(x, true); + d.reinit(x, true); + h.reinit(x, true); + + double gh,beta; + + // compute residual. if vector is + // zero, then short-circuit the + // full computation + if (!x.all_zero()) { - double eigen_beta_alpha = 0; + A.vmult(g,x); + g.add(-1.,b); + } + else + g.equ(-1.,b); + res = g.l2_norm(); - // resize the vectors, but do not set - // the values since they'd be overwritten - // soon anyway. - g.reinit(x, true); - d.reinit(x, true); - h.reinit(x, true); + conv = this->iteration_status(0, res, x); + if (conv != SolverControl::iterate) + return; - double gh,beta; + if (std::is_same::value == false) + { + preconditioner.vmult(h,g); - // compute residual. if vector is - // zero, then short-circuit the - // full computation - if (!x.all_zero()) - { - A.vmult(g,x); - g.add(-1.,b); - } - else - g.equ(-1.,b); - res = g.l2_norm(); + d.equ(-1.,h); + + gh = g*h; + } + else + { + d.equ(-1.,g); + gh = res*res; + } + + while (conv == SolverControl::iterate) + { + it++; + A.vmult(h,d); + + double alpha = d*h; + Assert(alpha != 0., ExcDivideByZero()); + alpha = gh/alpha; + + x.add(alpha,d); + res = std::sqrt(g.add_and_dot(alpha, h, g)); - conv = this->iteration_status(0, res, x); + print_vectors(it, x, g, d); + + conv = this->iteration_status(it, res, x); if (conv != SolverControl::iterate) - { - deallog.pop(); - return; - } + break; - if (std::is_same::value == false) + if (std::is_same::value + == false) { preconditioner.vmult(h,g); - d.equ(-1.,h); - - gh = g*h; + beta = gh; + Assert(beta != 0., ExcDivideByZero()); + gh = g*h; + beta = gh/beta; + d.sadd(beta,-1.,h); } else { - d.equ(-1.,g); + beta = gh; gh = res*res; + beta = gh/beta; + d.sadd(beta,-1.,g); } - while (conv == SolverControl::iterate) + this->coefficients_signal(alpha,beta); + // set up the vectors + // containing the diagonal + // and the off diagonal of + // the projected matrix. + if (do_eigenvalues) { - it++; - A.vmult(h,d); - - double alpha = d*h; - Assert(alpha != 0., ExcDivideByZero()); - alpha = gh/alpha; - - x.add(alpha,d); - res = std::sqrt(g.add_and_dot(alpha, h, g)); - - print_vectors(it, x, g, d); - - conv = this->iteration_status(it, res, x); - if (conv != SolverControl::iterate) - break; - - if (std::is_same::value - == false) - { - preconditioner.vmult(h,g); - - beta = gh; - Assert(beta != 0., ExcDivideByZero()); - gh = g*h; - beta = gh/beta; - d.sadd(beta,-1.,h); - } - else - { - beta = gh; - gh = res*res; - beta = gh/beta; - d.sadd(beta,-1.,g); - } - - this->coefficients_signal(alpha,beta); - // set up the vectors - // containing the diagonal - // and the off diagonal of - // the projected matrix. - if (do_eigenvalues) - { - diagonal.push_back(1./alpha + eigen_beta_alpha); - eigen_beta_alpha = beta/alpha; - offdiagonal.push_back(std::sqrt(beta)/alpha); - } - compute_eigs_and_cond(diagonal,offdiagonal,all_eigenvalues_signal, - all_condition_numbers_signal); + diagonal.push_back(1./alpha + eigen_beta_alpha); + eigen_beta_alpha = beta/alpha; + offdiagonal.push_back(std::sqrt(beta)/alpha); } + compute_eigs_and_cond(diagonal,offdiagonal,all_eigenvalues_signal, + all_condition_numbers_signal); } - catch (...) - { - deallog.pop(); - throw; - } + compute_eigs_and_cond(diagonal,offdiagonal,eigenvalues_signal, condition_number_signal); - deallog.pop(); - // in case of failure: throw exception if (conv != SolverControl::success) AssertThrow(false, SolverControl::NoConvergence (it, res)); diff --git a/include/deal.II/lac/solver_fire.h b/include/deal.II/lac/solver_fire.h index 316b664fb7..486f8dcce6 100644 --- a/include/deal.II/lac/solver_fire.h +++ b/include/deal.II/lac/solver_fire.h @@ -251,7 +251,7 @@ SolverFIRE::solve VectorType &x, const PreconditionerType &inverse_mass_matrix) { - deallog.push("FIRE"); + LogStream::Prefix prefix("FIRE"); // FIRE algorithm constants const double DELAYSTEP = 5; @@ -361,8 +361,6 @@ SolverFIRE::solve } // While we need to iterate. - deallog.pop(); - // In the case of failure: throw exception. if (conv != SolverControl::success) AssertThrow (false, diff --git a/include/deal.II/lac/solver_gmres.h b/include/deal.II/lac/solver_gmres.h index 4920926e28..ab8494d4df 100644 --- a/include/deal.II/lac/solver_gmres.h +++ b/include/deal.II/lac/solver_gmres.h @@ -777,7 +777,7 @@ SolverGMRES::solve (const MatrixType &A, //TODO:[?] Check, why there are two different start residuals. //TODO:[GK] Make sure the parameter in the constructor means maximum basis size - deallog.push("GMRES"); + LogStream::Prefix prefix("GMRES"); const unsigned int n_tmp_vectors = additional_data.max_n_tmp_vectors; // Generate an object where basis vectors are stored. @@ -1046,8 +1046,6 @@ SolverGMRES::solve (const MatrixType &A, if (!krylov_space_signal.empty()) krylov_space_signal(tmp_vectors); - deallog.pop(); - // in case of failure: throw exception AssertThrow(iteration_state == SolverControl::success, SolverControl::NoConvergence (accumulated_iterations, @@ -1172,7 +1170,7 @@ SolverFGMRES::solve (const MatrixType &A, const VectorType &b, const PreconditionerType &preconditioner) { - deallog.push("FGMRES"); + LogStream::Prefix prefix("FGMRES"); SolverControl::State iteration_state = SolverControl::iterate; @@ -1257,7 +1255,6 @@ SolverFGMRES::solve (const MatrixType &A, } while (iteration_state == SolverControl::iterate); - deallog.pop(); // in case of failure: throw exception if (iteration_state != SolverControl::success) AssertThrow(false, SolverControl::NoConvergence (accumulated_iterations, diff --git a/include/deal.II/lac/solver_minres.h b/include/deal.II/lac/solver_minres.h index d46cbb9d80..3c7f79566e 100644 --- a/include/deal.II/lac/solver_minres.h +++ b/include/deal.II/lac/solver_minres.h @@ -193,7 +193,7 @@ SolverMinRes::solve (const MatrixType &A, const VectorType &b, const PreconditionerType &preconditioner) { - deallog.push("minres"); + LogStream::Prefix prefix("minres"); // Memory allocation typename VectorMemory::Pointer Vu0 (this->memory); @@ -350,9 +350,6 @@ SolverMinRes::solve (const MatrixType &A, delta[1] = delta[2]; } - // Output - deallog.pop (); - // in case of failure: throw exception AssertThrow(conv == SolverControl::success, SolverControl::NoConvergence (j, r_l2)); diff --git a/include/deal.II/lac/solver_qmrs.h b/include/deal.II/lac/solver_qmrs.h index df6fef22f9..a3b40b1e7d 100644 --- a/include/deal.II/lac/solver_qmrs.h +++ b/include/deal.II/lac/solver_qmrs.h @@ -260,7 +260,7 @@ SolverQMRS::solve (const MatrixType &A, const VectorType &b, const PreconditionerType &preconditioner) { - deallog.push("QMRS"); + LogStream::Prefix prefix("QMRS"); // temporary vectors, allocated through the @p VectorMemory object at the // start of the actual solution process and deallocated at the end. @@ -291,9 +291,6 @@ SolverQMRS::solve (const MatrixType &A, } while (state.state == SolverControl::iterate); - // Output - deallog.pop(); - // in case of failure: throw exception AssertThrow(state.state == SolverControl::success, SolverControl::NoConvergence (step, diff --git a/include/deal.II/lac/solver_relaxation.h b/include/deal.II/lac/solver_relaxation.h index 0f212b37c8..7ae5880248 100644 --- a/include/deal.II/lac/solver_relaxation.h +++ b/include/deal.II/lac/solver_relaxation.h @@ -122,35 +122,26 @@ SolverRelaxation::solve (const MatrixType &A, VectorType &d = *Vd; d.reinit(x); - deallog.push("Relaxation"); + LogStream::Prefix prefix("Relaxation"); int iter=0; - try + // Main loop + for (; conv==SolverControl::iterate; iter++) { - // Main loop - for (; conv==SolverControl::iterate; iter++) - { - // Compute residual - A.vmult(r,x); - r.sadd(-1.,1.,b); - - // The required norm of the - // (preconditioned) - // residual is computed in - // criterion() and stored - // in res. - conv = this->iteration_status (iter, r.l2_norm(), x); - if (conv != SolverControl::iterate) - break; - R.step(x,b); - } + // Compute residual + A.vmult(r,x); + r.sadd(-1.,1.,b); + + // The required norm of the + // (preconditioned) + // residual is computed in + // criterion() and stored + // in res. + conv = this->iteration_status (iter, r.l2_norm(), x); + if (conv != SolverControl::iterate) + break; + R.step(x,b); } - catch (...) - { - deallog.pop(); - throw; - } - deallog.pop(); // in case of failure: throw exception AssertThrow(conv == SolverControl::success, diff --git a/include/deal.II/lac/solver_richardson.h b/include/deal.II/lac/solver_richardson.h index 487f8e5590..23e7b7fb9c 100644 --- a/include/deal.II/lac/solver_richardson.h +++ b/include/deal.II/lac/solver_richardson.h @@ -217,40 +217,30 @@ SolverRichardson::solve (const MatrixType &A, VectorType &d = *Vd; d.reinit(x); - deallog.push("Richardson"); + LogStream::Prefix prefix("Richardson"); - try + // Main loop + while (conv==SolverControl::iterate) { - // Main loop - while (conv==SolverControl::iterate) - { - // Do not use residual, - // but do it in 2 steps - A.vmult(r,x); - r.sadd(-1.,1.,b); - preconditioner.vmult(d,r); - - // get the required norm of the (possibly preconditioned) - // residual - last_criterion = criterion(r, d); - conv = this->iteration_status (iter, last_criterion, x); - if (conv != SolverControl::iterate) - break; - - x.add(additional_data.omega,d); - print_vectors(iter,x,r,d); - - ++iter; - } - } - catch (...) - { - deallog.pop(); - throw; + // Do not use residual, + // but do it in 2 steps + A.vmult(r,x); + r.sadd(-1.,1.,b); + preconditioner.vmult(d,r); + + // get the required norm of the (possibly preconditioned) + // residual + last_criterion = criterion(r, d); + conv = this->iteration_status (iter, last_criterion, x); + if (conv != SolverControl::iterate) + break; + + x.add(additional_data.omega,d); + print_vectors(iter,x,r,d); + + ++iter; } - deallog.pop(); - // in case of failure: throw exception if (conv != SolverControl::success) AssertThrow(false, SolverControl::NoConvergence (iter, @@ -284,37 +274,28 @@ SolverRichardson::Tsolve (const MatrixType &A, VectorType &d = *Vd; d.reinit(x); - deallog.push("RichardsonT"); + LogStream::Prefix prefix("RichardsonT"); - try - { - // Main loop - while (conv==SolverControl::iterate) - { - // Do not use Tresidual, - // but do it in 2 steps - A.Tvmult(r,x); - r.sadd(-1.,1.,b); - preconditioner.Tvmult(d,r); - - last_criterion = criterion(r, d); - conv = this->iteration_status (iter, last_criterion, x); - if (conv != SolverControl::iterate) - break; - - x.add(additional_data.omega,d); - print_vectors(iter,x,r,d); - - ++iter; - } - } - catch (...) + // Main loop + while (conv==SolverControl::iterate) { - deallog.pop(); - throw; + // Do not use Tresidual, + // but do it in 2 steps + A.Tvmult(r,x); + r.sadd(-1.,1.,b); + preconditioner.Tvmult(d,r); + + last_criterion = criterion(r, d); + conv = this->iteration_status (iter, last_criterion, x); + if (conv != SolverControl::iterate) + break; + + x.add(additional_data.omega,d); + print_vectors(iter,x,r,d); + + ++iter; } - deallog.pop(); // in case of failure: throw exception if (conv != SolverControl::success) AssertThrow(false, SolverControl::NoConvergence (iter, last_criterion));