From 1b3834d7b59da2e809d320992c632697355d63b6 Mon Sep 17 00:00:00 2001 From: Pauli Virtanen Date: Fri, 12 Apr 2013 18:16:31 +0300 Subject: BUG: linalg: make umath_linalg to track errors from all inner loop iterations This ensures that the FP invalid flag always reflects the return code from LAPACK. Fixes a bug in 63a8aef81 where umath_linalg raises a warning only if the error occurs in the last iteration of the ufunc inner loop. --- numpy/linalg/umath_linalg.c.src | 85 +++++++++++++++++------------------------ 1 file changed, 34 insertions(+), 51 deletions(-) (limited to 'numpy/linalg/umath_linalg.c.src') diff --git a/numpy/linalg/umath_linalg.c.src b/numpy/linalg/umath_linalg.c.src index c9d55f635..eadbde8e7 100644 --- a/numpy/linalg/umath_linalg.c.src +++ b/numpy/linalg/umath_linalg.c.src @@ -374,10 +374,23 @@ offset_ptr(void* ptr, ptrdiff_t offset) return (void*)((npy_uint8*)ptr + offset); } +static inline int +get_fp_invalid_and_clear() +{ + int status; + status = PyUFunc_getfperr(); + return !!(status & UFUNC_FPE_INVALID); +} + static inline void -clear_fp_errors() +set_fp_invalid_or_clear(int error_occurred) { - PyUFunc_getfperr(); + if (error_occurred) { + npy_set_floatstatus_invalid(); + } + else { + PyUFunc_getfperr(); + } } /* @@ -862,10 +875,6 @@ nan_@TYPE@_matrix(void *dst_in, const LINEARIZE_DATA_t* data) } dst += data->row_strides/sizeof(@typ@); } - - /* if this function gets called, the module generates nans */ - /* set the fpe status to invalid to notify that nans were generated */ - npy_set_floatstatus_invalid(); } static inline void @@ -1786,7 +1795,7 @@ static inline void size_t outer_dim = *dimensions++; size_t op_count = (JOBZ=='N')?2:3; EIGH_PARAMS_t eigh_params; - int error_occurred = 0; + int error_occurred = get_fp_invalid_and_clear(); for (iter=0; iter < op_count; ++iter) { outer_steps[iter] = (ptrdiff_t) steps[iter]; @@ -1843,10 +1852,7 @@ static inline void release_@lapack_func@(&eigh_params); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } /**end repeat**/ @@ -1973,7 +1979,7 @@ static void { GESV_PARAMS_t params; fortran_int n, nrhs; - int error_occurred = 0; + int error_occurred = get_fp_invalid_and_clear(); INIT_OUTER_LOOP_3 n = (fortran_int)dimensions[0]; @@ -2001,10 +2007,7 @@ static void release_@lapack_func@(¶ms); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } static void @@ -2012,7 +2015,7 @@ static void void *NPY_UNUSED(func)) { GESV_PARAMS_t params; - int error_occurred; + int error_occurred = get_fp_invalid_and_clear(); fortran_int n; INIT_OUTER_LOOP_3 @@ -2039,10 +2042,7 @@ static void release_@lapack_func@(¶ms); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } static void @@ -2051,7 +2051,7 @@ static void { GESV_PARAMS_t params; fortran_int n; - int error_occurred = 0; + int error_occurred = get_fp_invalid_and_clear(); INIT_OUTER_LOOP_2 n = (fortran_int)dimensions[0]; @@ -2076,11 +2076,9 @@ static void release_@lapack_func@(¶ms); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } + /**end repeat**/ @@ -2150,7 +2148,7 @@ static void @TYPE@_cholesky(char uplo, char **args, npy_intp *dimensions, npy_intp *steps) { POTR_PARAMS_t params; - int error_occurred = 0; + int error_occurred = get_fp_invalid_and_clear(); fortran_int n; INIT_OUTER_LOOP_2 @@ -2175,10 +2173,7 @@ static void release_@lapack_func@(¶ms); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } static void @@ -2592,7 +2587,7 @@ static inline void size_t iter; size_t outer_dim = *dimensions++; size_t op_count = 2; - int error_occurred = 0; + int error_occurred = get_fp_invalid_and_clear(); GEEV_PARAMS_t geev_params; STACK_TRACE; @@ -2668,10 +2663,7 @@ static inline void release_@lapack_func@(&geev_params); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } static void @@ -3038,7 +3030,7 @@ static inline void npy_intp* steps) { ptrdiff_t outer_steps[3]; - int error_occurred = 0; + int error_occurred = get_fp_invalid_and_clear(); size_t iter; size_t outer_dim = *dimensions++; size_t op_count = (JOBZ=='N')?2:4; @@ -3110,10 +3102,7 @@ static inline void release_@lapack_func@(¶ms); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } /**end repeat*/ @@ -3519,7 +3508,7 @@ static void ) { POTRS_PARAMS_t params; - int error_occurred = 0; + int error_occurred = get_fp_invalid_and_clear(); fortran_int n, nrhs; INIT_OUTER_LOOP_3 @@ -3553,10 +3542,7 @@ static void release_@lapack_func@(¶ms); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } static void @@ -3679,7 +3665,7 @@ static inline void npy_intp *steps) { POTRI_PARAMS_t params; - int error_occurred = 0; + int error_occurred = get_fp_invalid_and_clear(); fortran_int n; INIT_OUTER_LOOP_2 @@ -3706,10 +3692,7 @@ static inline void release_@potri@(¶ms); } - if (!error_occurred) { - /* lapack success, clear FP error status */ - clear_fp_errors(); - } + set_fp_invalid_or_clear(error_occurred); } static void -- cgit v1.2.1