Mercurial > octave
changeset 33203:371217bc9733
Replace obscure Fortran coding pattern in liboctave with modern C++ (bug #65343)
* CMatrix.cc, CSparse.cc, dMatrix.cc, dSparse.cc, fCMatrix.cc, fMatrix.cc:
Add #include for <cmath> and <limits>. Replace all instances of
"rcond_plus_one == 1.0" with "is_singular (rcond)". Delete variable
instantiation for "rcond_plus_one".
* CMatrix.cc, CSparse.cc, dMatrix.cc, dSparse.cc, fCMatrix.cc, fMatrix.cc
(is_singular): New static inline bool function to check whether a matrix
is singular by checking the size of the reciprocal condition number.
| author | Rik <rik@octave.org> |
|---|---|
| date | Wed, 13 Mar 2024 14:11:56 -0700 |
| parents | cab636a9c37b |
| children | ccd0c0eaa323 |
| files | liboctave/array/CMatrix.cc liboctave/array/CSparse.cc liboctave/array/dMatrix.cc liboctave/array/dSparse.cc liboctave/array/fCMatrix.cc liboctave/array/fMatrix.cc |
| diffstat | 6 files changed, 115 insertions(+), 236 deletions(-) [+] |
line wrap: on
line diff
--- a/liboctave/array/CMatrix.cc Tue Mar 12 22:07:22 2024 -0700 +++ b/liboctave/array/CMatrix.cc Wed Mar 13 14:11:56 2024 -0700 @@ -28,6 +28,7 @@ #endif #include <algorithm> +#include <cmath> #include <complex> #include <istream> #include <limits> @@ -733,6 +734,14 @@ return anorm; } +// Local function to check if matrix is singular based on rcond. +static inline +bool +is_singular (const double rcond) +{ + return (std::abs (rcond) <= std::numeric_limits<double>::epsilon ()); +} + ComplexMatrix ComplexMatrix::inverse () const { @@ -1603,10 +1612,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1702,10 +1708,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1795,10 +1798,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1891,10 +1891,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { if (sing_handler) sing_handler (rcon);
--- a/liboctave/array/CSparse.cc Tue Mar 12 22:07:22 2024 -0700 +++ b/liboctave/array/CSparse.cc Wed Mar 13 14:11:56 2024 -0700 @@ -27,8 +27,10 @@ # include "config.h" #endif +#include <cmath> #include <complex> #include <istream> +#include <limits> #include <ostream> #include "quit.h" @@ -651,6 +653,14 @@ return retval; } +// Local function to check if matrix is singular based on rcond. +static inline +bool +is_singular (const double rcond) +{ + return (std::abs (rcond) <= std::numeric_limits<double>::epsilon ()); +} + SparseComplexMatrix SparseComplexMatrix::inverse () const { @@ -1723,10 +1733,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -2006,10 +2013,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -2236,10 +2240,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -2519,10 +2520,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -2769,10 +2767,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -3071,10 +3066,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -3324,10 +3316,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -3625,10 +3614,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4413,10 +4399,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4547,10 +4530,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4684,10 +4664,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4851,10 +4828,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5027,10 +5001,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5158,10 +5129,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5297,10 +5265,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5468,10 +5433,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5633,11 +5595,8 @@ else rcond = 1.; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - if (status == UMFPACK_WARNING_singular_matrix - || rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + || is_singular (rcond) || octave::math::isnan (rcond)) { UMFPACK_ZNAME (report_numeric) (Numeric, control); @@ -5784,10 +5743,7 @@ } else { - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -6034,10 +5990,7 @@ } else { - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -6317,10 +6270,7 @@ } else { - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -6546,10 +6496,7 @@ } else { - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -6682,11 +6629,8 @@ retval.maybe_compress (); rcond = Info (UMFPACK_RCOND); - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - if (status == UMFPACK_WARNING_singular_matrix - || rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + || is_singular (rcond) || octave::math::isnan (rcond)) { err = -2;
--- a/liboctave/array/dMatrix.cc Tue Mar 12 22:07:22 2024 -0700 +++ b/liboctave/array/dMatrix.cc Wed Mar 13 14:11:56 2024 -0700 @@ -28,6 +28,7 @@ #endif #include <algorithm> +#include <cmath> #include <istream> #include <limits> #include <ostream> @@ -447,6 +448,14 @@ return anorm; } +// Local function to check if matrix is singular based on rcond. +static inline +bool +is_singular (const double rcond) +{ + return (std::abs (rcond) <= std::numeric_limits<double>::epsilon ()); +} + Matrix Matrix::inverse () const { @@ -1274,10 +1283,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1372,10 +1378,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1461,10 +1464,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1551,10 +1551,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { if (sing_handler) sing_handler (rcon);
--- a/liboctave/array/dSparse.cc Tue Mar 12 22:07:22 2024 -0700 +++ b/liboctave/array/dSparse.cc Wed Mar 13 14:11:56 2024 -0700 @@ -27,7 +27,9 @@ # include "config.h" #endif +#include <cmath> #include <istream> +#include <limits> #include <ostream> #include "quit.h" @@ -594,6 +596,14 @@ */ +// Local function to check if matrix is singular based on rcond. +static inline +bool +is_singular (const double rcond) +{ + return (std::abs (rcond) <= std::numeric_limits<double>::epsilon ()); +} + SparseMatrix SparseMatrix::inverse () const { @@ -1657,10 +1667,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -1940,10 +1947,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -2172,10 +2176,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -2457,10 +2458,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -2711,10 +2709,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -3012,10 +3007,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -3267,10 +3259,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -3570,10 +3559,7 @@ octave::warn_singular_matrix (rcond); } - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4363,10 +4349,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4496,10 +4479,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4633,10 +4613,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4800,10 +4777,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -4977,10 +4951,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5142,10 +5113,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5304,10 +5272,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5498,10 +5463,7 @@ if (err != 0) err = -2; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -5674,11 +5636,8 @@ else rcond = 1.; - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - if (status == UMFPACK_WARNING_singular_matrix - || rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + || is_singular (rcond) || octave::math::isnan (rcond)) { UMFPACK_DNAME (report_numeric) (Numeric, control); @@ -5826,10 +5785,7 @@ } else { - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -6048,10 +6004,7 @@ } else { - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -6302,10 +6255,7 @@ } else { - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2; @@ -6545,10 +6495,7 @@ } else { - // Prevent use of extra precision. - double rcond_plus_one = rcond + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcond)) + if (is_singular (rcond) || octave::math::isnan (rcond)) { err = -2;
--- a/liboctave/array/fCMatrix.cc Tue Mar 12 22:07:22 2024 -0700 +++ b/liboctave/array/fCMatrix.cc Wed Mar 13 14:11:56 2024 -0700 @@ -28,6 +28,7 @@ #endif #include <algorithm> +#include <cmath> #include <complex> #include <istream> #include <limits> @@ -736,6 +737,14 @@ return anorm; } +// Local function to check if matrix is singular based on rcond. +static inline +bool +is_singular (const float rcond) +{ + return (std::abs (rcond) <= std::numeric_limits<float>::epsilon ()); +} + FloatComplexMatrix FloatComplexMatrix::inverse () const { @@ -1605,10 +1614,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - float rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1708,10 +1714,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - float rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1805,10 +1808,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - float rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1901,10 +1901,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - float rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { if (sing_handler) sing_handler (rcon);
--- a/liboctave/array/fMatrix.cc Tue Mar 12 22:07:22 2024 -0700 +++ b/liboctave/array/fMatrix.cc Wed Mar 13 14:11:56 2024 -0700 @@ -28,6 +28,7 @@ #endif #include <algorithm> +#include <cmath> #include <istream> #include <limits> #include <ostream> @@ -453,6 +454,14 @@ return anorm; } +// Local function to check if matrix is singular based on rcond. +static inline +bool +is_singular (const float rcond) +{ + return (std::abs (rcond) <= std::numeric_limits<float>::epsilon ()); +} + FloatMatrix FloatMatrix::inverse () const { @@ -1284,10 +1293,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - float rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1388,10 +1394,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - float rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1485,10 +1488,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - float rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { info = -2; @@ -1572,10 +1572,7 @@ if (info != 0) info = -2; - // Prevent use of extra precision. - float rcond_plus_one = rcon + 1.0; - - if (rcond_plus_one == 1.0 || octave::math::isnan (rcon)) + if (is_singular (rcon) || octave::math::isnan (rcon)) { if (sing_handler) sing_handler (rcon);