Avoid calling *trsv and *gemv with invalid parameter values#174
Avoid calling *trsv and *gemv with invalid parameter values#174rgommers wants to merge 4 commits into
*trsv and *gemv with invalid parameter values#174Conversation
This upstreams a fix that SciPy has carried for a long time. Note that in `dgstrf.c`, the call to `dsnode_bmod` comes before the call to `dpivotL` where a singularity is detected, so the fix can not be lifted up to that file.
When factorizing very sparse or singular matrices, the leading dimension (nsupr) can be smaller than the segment size (segsze), which results in invalid parameters being passed to `trsv`. Add a guard to abort early with a clear error message instead. Upstreamed from SciPy, which has carried this patch for a long time.
The ILU pivot routines had two singularity checks with commented-out fprintf/exit calls that silently continued. Replace these with ABORT() to fail with a clear error message instead. Upstreamed from SciPy, which has carried this patch for a long time.
…ingular matrix When the matrix is singular (pivmax == 0.0), pivptr may be out of bounds for lsub_ptr (pivptr >= nsupr). Add a guard to use diagind as the pivot row when pivptr is out of range. Upstreamed from SciPy, which has carried this patch for a long time. Also cleans up some code that was ifdef'd out already.
|
I like the changes! A potential further improvement would be to extend the message "failed to factorize matrix" by reason given in the commit desciptions. @xiaoyeli What do you think? |
|
Gentle ping @xiaoyeli , would you have time to take a look at this PR? |
|
Gentle ping @xiaoyeli , could you have a look? |
|
Thanks for looking at this @xiaoyeli. I don't think the performance impact is measurable at all, and certainly it cannot be "a lot". It's a single integer comparison (one cycle at most, <1 ns) for both changes you point at, followed by an O(N^2) TRSV call (and in the first patch, TRSV + GEMV), and it's a not-taken branch, so the runtime cost is noise. Conceptually though, I agree it'd be nice to lift up robustness checks as far as possible. However, I spent significant time on it after your suggestion to do that, and wasn't able to come up with a clean way of doing it, as commented on in #167 (comment). At a higher level, even if succcessful, it would increase the diff size I think, because the |
|
It is not just the comparison, it is the conditional jump that might invalidate the whole pipeline. I am not aware of a way to hint the compiler in C with the likelihood of a branch (unlike C++20 with Can you run some of the tests and see if your patch has an impact on the performance? |
|
A conditional jump: true in principle indeed. I don't believe it matters here because (a) it's a forward branch, and CPUs default to not-taken prediction when cold, (b) the condition is essentially always false (modulo the crash we're trying to avoid), which is as easy as it gets for dynamic prediction - no modern CPU should get this wrong, and (c) even if they get it wrong once, we're then taking about a cost that's maybe 20 cycles, still far below the cost of a single
You can define a macro in C built on That shouldn't affect branch prediction on modern compilers, but it helps by ensuring the the generated code for the likely path is inline. It seems unlikely that it's going to do have any effect here, but I'd be willing to add it, and perhaps you see uses for it elsewhere? This would be the easy version: #if defined(__GNUC__) /* covers GCC, Clang, ICX, and most GCC-compatible compilers */
# define EXPECT_UNLIKELY(x) __builtin_expect(!!(x), 0)
#else
# define EXPECT_UNLIKELY(x) (x)
#endif
#endifAn actual compile check in CMake/Make is the slightly more thorough version if you care about niche compilers.
I can do that - it will take me some time though. It's almost certainly going to yield overlapping error bars though, even when sampling a lot. I suspect a check with |
This upstreams patches that SciPy has carried for a long time, as discussed in gh-167. There are four separate fixes - this PR has them all, one commit per issue. The commit messages have more detailed explanation.
I also investigated moving the fixes around as suggested in gh-167, but that doesn't quite work as commented on in #167 (comment).