From: kanschat Date: Fri, 15 Oct 2010 14:09:00 +0000 (+0000) Subject: catch more divide by zero and infinite numbers X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ac1610982f847fd95e2c1858e3dc27f0b34fad4f;p=dealii-svn.git catch more divide by zero and infinite numbers git-svn-id: https://svn.dealii.org/trunk@22349 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/source/flow_function.cc b/deal.II/base/source/flow_function.cc index 535b4b916c..12cf8a8209 100644 --- a/deal.II/base/source/flow_function.cc +++ b/deal.II/base/source/flow_function.cc @@ -185,7 +185,9 @@ namespace Functions const double Re) : radius(r), Reynolds(Re) - {} + { + Assert(Reynolds != 0., ExcMessage("Reynolds number cannot be zero")); + } template diff --git a/deal.II/lac/include/lac/filtered_matrix.h b/deal.II/lac/include/lac/filtered_matrix.h index 20de47e91c..2bc9a6697f 100644 --- a/deal.II/lac/include/lac/filtered_matrix.h +++ b/deal.II/lac/include/lac/filtered_matrix.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors +// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -922,22 +922,27 @@ FilteredMatrix::apply_constraints ( const bool /* matrix_is_symmetric */) const { GrowingVectorMemory mem; - VECTOR* tmp_vector = mem.alloc(); + typename VectorMemory::Pointer tmp_vector(mem); tmp_vector->reinit(v); const_index_value_iterator i = constraints.begin(); const const_index_value_iterator e = constraints.end(); for (; i!=e; ++i) - (*tmp_vector)(i->first) = -i->second; + { + Assert(numbers::is_finite(i->second), ExcNumberNotFinite()); + (*tmp_vector)(i->first) = -i->second; + } // This vmult is without bc, to get // the rhs correction in a correct // way. matrix->vmult_add(v, *tmp_vector); - mem.free(tmp_vector); // finally set constrained // entries themselves for (i=constraints.begin(); i!=e; ++i) - v(i->first) = i->second; + { + Assert(numbers::is_finite(i->second), ExcNumberNotFinite()); + v(i->first) = i->second; + } } @@ -968,7 +973,10 @@ FilteredMatrix::post_filter (const VECTOR &in, const_index_value_iterator i = constraints.begin(); const const_index_value_iterator e = constraints.end(); for (; i!=e; ++i) - out(i->first) = in(i->first); + { + Assert(numbers::is_finite(in(i->first)), ExcNumberNotFinite()); + out(i->first) = in(i->first); + } } diff --git a/deal.II/lac/include/lac/sparse_matrix.templates.h b/deal.II/lac/include/lac/sparse_matrix.templates.h index 9d3f561727..f107bdc6c9 100644 --- a/deal.II/lac/include/lac/sparse_matrix.templates.h +++ b/deal.II/lac/include/lac/sparse_matrix.templates.h @@ -1342,6 +1342,7 @@ SparseMatrix::precondition_SSOR (Vector &dst, // divide by diagonal element *dst_ptr -= s * om; + Assert(val[*rowstart_ptr]!= 0., ExcDivideByZero()); *dst_ptr /= val[*rowstart_ptr]; }; @@ -1363,6 +1364,7 @@ SparseMatrix::precondition_SSOR (Vector &dst, s += val[j] * dst(cols->colnums[j]); *dst_ptr -= s * om; + Assert(val[*rowstart_ptr]!= 0., ExcDivideByZero()); *dst_ptr /= val[*rowstart_ptr]; }; return; @@ -1395,6 +1397,7 @@ SparseMatrix::precondition_SSOR (Vector &dst, // divide by diagonal element *dst_ptr -= s * om; + Assert(val[*rowstart_ptr]!= 0., ExcDivideByZero()); *dst_ptr /= val[*rowstart_ptr]; }; @@ -1418,6 +1421,7 @@ SparseMatrix::precondition_SSOR (Vector &dst, for (unsigned int j=first_right_of_diagonal_index; jcolnums[j]); *dst_ptr -= s * om; + Assert(val[*rowstart_ptr]!= 0., ExcDivideByZero()); *dst_ptr /= val[*rowstart_ptr]; }; } @@ -1481,7 +1485,8 @@ SparseMatrix::SOR (Vector& dst, if (col < row) s -= val[j] * dst(col); } - + + Assert(val[cols->rowstart[row]]!= 0., ExcDivideByZero()); dst(row) = s * om / val[cols->rowstart[row]]; } } @@ -1508,6 +1513,7 @@ SparseMatrix::TSOR (Vector& dst, if (cols->colnums[j] > row) s -= val[j] * dst(cols->colnums[j]); + Assert(val[cols->rowstart[row]]!= 0., ExcDivideByZero()); dst(row) = s * om / val[cols->rowstart[row]]; if (row == 0) @@ -1551,6 +1557,7 @@ SparseMatrix::PSOR (Vector& dst, } } + Assert(val[cols->rowstart[row]]!= 0., ExcDivideByZero()); dst(row) = s * om / val[cols->rowstart[row]]; } } @@ -1587,6 +1594,7 @@ SparseMatrix::TPSOR (Vector& dst, s -= val[j] * dst(col); } + Assert(val[cols->rowstart[row]]!= 0., ExcDivideByZero()); dst(row) = s * om / val[cols->rowstart[row]]; } } @@ -1647,6 +1655,7 @@ SparseMatrix::SOR_step (Vector &v, { s -= val[j] * v(cols->colnums[j]); } + Assert(val[cols->rowstart[row]]!= 0., ExcDivideByZero()); v(row) += s * om / val[cols->rowstart[row]]; } } @@ -1675,6 +1684,7 @@ SparseMatrix::TSOR_step (Vector &v, { s -= val[j] * v(cols->colnums[j]); } + Assert(val[cols->rowstart[row]]!= 0., ExcDivideByZero()); v(row) += s * om / val[cols->rowstart[row]]; } } @@ -1726,6 +1736,7 @@ SparseMatrix::SSOR (Vector& dst, } } dst(i) -= s * om; + Assert(val[cols->rowstart[i]]!= 0., ExcDivideByZero()); dst(i) /= val[cols->rowstart[i]]; } @@ -1740,6 +1751,7 @@ SparseMatrix::SSOR (Vector& dst, if (static_cast(i)rowstart[i]]!= 0., ExcDivideByZero()); dst(i) -= s * om / val[cols->rowstart[i]]; } }