From: kanschat Date: Thu, 21 Oct 2010 14:08:18 +0000 (+0000) Subject: more tests on infinite and divide by zero X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=656c1aa3589f55f3b67ee3d48942735f4fec7fa7;p=dealii-svn.git more tests on infinite and divide by zero git-svn-id: https://svn.dealii.org/trunk@22416 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/solver_cg.h b/deal.II/lac/include/lac/solver_cg.h index 3f182e64da..7a561e9e6f 100644 --- a/deal.II/lac/include/lac/solver_cg.h +++ b/deal.II/lac/include/lac/solver_cg.h @@ -367,6 +367,7 @@ SolverCG::solve (const MATRIX &A, A.vmult(Ad,d); alpha = d*Ad; + Assert(alpha != 0., ExcDivideByZero()); alpha = gh/alpha; g.add(alpha,Ad); @@ -382,6 +383,7 @@ SolverCG::solve (const MATRIX &A, precondition.vmult(h,g); beta = gh; + Assert(beta != 0., ExcDivideByZero()); gh = g*h; beta = gh/beta; @@ -397,7 +399,7 @@ SolverCG::solve (const MATRIX &A, eigen_beta_alpha = beta/alpha; offdiagonal.push_back(std::sqrt(beta)/alpha); } - + if (additional_data.compute_all_condition_numbers && (diagonal.size()>1)) { TridiagonalMatrix T(diagonal.size(), true); diff --git a/deal.II/lac/include/lac/vector.templates.h b/deal.II/lac/include/lac/vector.templates.h index 80819f7d55..858b78cdbd 100644 --- a/deal.II/lac/include/lac/vector.templates.h +++ b/deal.II/lac/include/lac/vector.templates.h @@ -393,6 +393,8 @@ Number Vector::operator * (const Vector& v) const sum3 += *X++ * Number(numbers::NumberTraits::conjugate(*Y++)); sum += sum3; + Assert(numbers::is_finite(sum), ExcNumberNotFinite()); + return sum; } @@ -452,6 +454,8 @@ Vector::norm_sqr () const sum3 += numbers::NumberTraits::abs_square(*X++); sum += sum3; + Assert(numbers::is_finite(sum), ExcNumberNotFinite()); + return sum; } @@ -855,8 +859,7 @@ template void Vector::equ (const Number a, const Vector& u) { - Assert (numbers::is_finite(a), - ExcMessage("The given value is not finite but either infinite or Not A Number (NaN)")); + Assert (numbers::is_finite(a), ExcNumberNotFinite()); Assert (vec_size!=0, ExcEmptyObject()); Assert (vec_size == u.vec_size, ExcDimensionMismatch(vec_size, u.vec_size)); @@ -879,8 +882,7 @@ template void Vector::equ (const Number a, const Vector& u) { - Assert (numbers::is_finite(a), - ExcMessage("The given value is not finite but either infinite or Not A Number (NaN)")); + Assert (numbers::is_finite(a), ExcNumberNotFinite()); Assert (vec_size!=0, ExcEmptyObject()); Assert (vec_size == u.vec_size, ExcDimensionMismatch(vec_size, u.vec_size));