From 31dbde380785be32893bfd3bba03df08c5598439 Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Wed, 15 Jun 2016 11:14:49 +0200 Subject: [PATCH] parpack: switch solver related Asserts to AssertThrow; minor clean-up of in-code comments --- include/deal.II/lac/parpack_solver.h | 43 ++++++++++++++-------------- 1 file changed, 22 insertions(+), 21 deletions(-) diff --git a/include/deal.II/lac/parpack_solver.h b/include/deal.II/lac/parpack_solver.h index 5cfbf42591..2e5f0058b6 100644 --- a/include/deal.II/lac/parpack_solver.h +++ b/include/deal.II/lac/parpack_solver.h @@ -543,7 +543,6 @@ internal_reinit(const IndexSet &locally_owned_dofs) ldv = nloc; v.resize (ldv*ncv, 0.0); - // TODO: add optional input for resid resid.resize(nloc, 1.0); // work arrays for ARPACK @@ -715,6 +714,7 @@ void PArpackSolver::solve int nev = n_eigenvalues; int n_inside_arpack = nloc; + // IDO = 99: done while (ido != 99) { // call of ARPACK pdnaupd routine @@ -737,10 +737,10 @@ void PArpackSolver::solve { switch (ido) { -// compute Y = OP * X where -// IPNTR(1) is the pointer into WORKD for X, -// IPNTR(2) is the pointer into WORKD for Y. case -1: + // compute Y = OP * X where + // IPNTR(1) is the pointer into WORKD for X, + // IPNTR(2) is the pointer into WORKD for Y. { const int shift_x = ipntr[0]-1; const int shift_y = ipntr[1]-1; @@ -767,13 +767,13 @@ void PArpackSolver::solve } break; -// compute Y = OP * X where -// IPNTR(1) is the pointer into WORKD for X, -// IPNTR(2) is the pointer into WORKD for Y. -// In mode 3,4 and 5, the vector B * X is already -// available in WORKD(ipntr(3)). It does not -// need to be recomputed in forming OP * X. case 1: + // compute Y = OP * X where + // IPNTR(1) is the pointer into WORKD for X, + // IPNTR(2) is the pointer into WORKD for Y. + // In mode 3,4 and 5, the vector B * X is already + // available in WORKD(ipntr(3)). It does not + // need to be recomputed in forming OP * X. { const int shift_x = ipntr[0]-1; const int shift_y = ipntr[1]-1; @@ -812,10 +812,10 @@ void PArpackSolver::solve } break; -// compute Y = B * X where -// IPNTR(1) is the pointer into WORKD for X, -// IPNTR(2) is the pointer into WORKD for Y. case 2: + // compute Y = B * X where + // IPNTR(1) is the pointer into WORKD for X, + // IPNTR(2) is the pointer into WORKD for Y. { const int shift_x = ipntr[0]-1; @@ -843,20 +843,20 @@ void PArpackSolver::solve break; default: - Assert (false, PArpackExcIdo(ido)); + AssertThrow (false, PArpackExcIdo(ido)); break; } } break; default: - Assert (false, PArpackExcMode(mode)); + AssertThrow (false, PArpackExcMode(mode)); break; } } if (info<0) { - Assert (false, PArpackExcInfoPdnaupd(info)); + AssertThrow (false, PArpackExcInfoPdnaupd(info)); } else { @@ -889,15 +889,15 @@ void PArpackSolver::solve if (info == 1) { - Assert (false, PArpackExcInfoMaxIt(control().max_steps())); + AssertThrow (false, PArpackExcInfoMaxIt(control().max_steps())); } else if (info == 3) { - Assert (false, PArpackExcNoShifts(1)); + AssertThrow (false, PArpackExcNoShifts(1)); } else if (info!=0) { - Assert (false, PArpackExcInfoPdneupd(info)); + AssertThrow (false, PArpackExcInfoPdneupd(info)); } for (size_type i=0; i::solve eigenvalues_im[i]); } - Assert (iparam[4] == (int)n_eigenvalues, - PArpackExcConvergedEigenvectors(n_eigenvalues,iparam[4])); + // Throw an error if the solver did not converge. + AssertThrow (iparam[4] == (int)n_eigenvalues, + PArpackExcConvergedEigenvectors(n_eigenvalues,iparam[4])); // both PDNAUPD and PDSAUPD compute eigenpairs of inv[A - sigma*M]*M // with respect to a semi-inner product defined by M. -- 2.39.5