Improve support for linear algebra

[OlivierHnt]

Some linear algebra functions are accessible from Base, e.g. * and inv.
This PR implements methods for these functions as part of the LinearAlgebra extension to avoid getting error messages and the need to load the more specialised package IntervalLinearAlgebra.

Specifically,

• adds a simple matrix inversion algorithm based on the contraction mapping theorem (x-ref: [enhancement]: inv and det (needed by IntervalRootFinding) IntervalLinearAlgebra.jl#134)
• implements Rump's algorithm (probably the most naive one, which does not rely on OpenBlas to avoid multithreading issues, x-ref: [enhancement]: Ship a correctly rounded threaded OpenBLAS as an Artifact IntervalLinearAlgebra.jl#131)

It also defines some methods for the operator norm (opnorm) since the default algorithm used the syntax sum::Tsum = 0 which leads to an unnecessary "NG" flag.

I took the tests for matrix multiplication from IntervalLinearAlgebra.

cc @lucaferranti @orkolorko

[schillic]

The matrix multiplication has become a lot more conservative. Was this intended?

julia> A = [interval(2, 4) interval(-2, 1); interval(-1, 2) interval(2, 4)];

julia> A*A  # v0.22.16
2×2 Matrix{Interval{Float64}}:
  [0.0, 18.0]_com  [-16.0, 8.0]_com
 [-8.0, 16.0]_com    [0.0, 18.0]_com

julia> A*A  # v0.22.17
2×2 Matrix{Interval{Float64}}:
  [-2.0, 19.5001]_com  [-16.0, 10.0]_com
 [-10.0, 16.0]_com      [-2.0, 19.5001]_com

I think this change may have affected the usefulness of the package quite a bit.
Is there a way to obtain the previous results?
If so, is there a way to change the behavior globally?
If not, is that planned?

[OlivierHnt]

Yes this is intended, the performance is really bad otherwise.
We could improve the current situation by adding a matrix_mode function controlling whether to use :fast (default) and slow, similar to IntervalArithmetic.power_mode and PowerMode for powers.

[schillic]

We could improve the current situation by adding a matrix_mode function controlling whether to use :fast (default) and slow, similar to IntervalArithmetic.power_mode and PowerMode for powers.

I think that would be good.
If you can do a computation faster - great.
If you can get a tighter (but still sound) result - great.
But if you have to make a trade-off, I would say that the user should have the choice. Some users expect the textbook definition of interval arithmetic, other users want tight results, other users want fast results. Of course you can make a default choice.

(I also think that the latest release should have been marked as a bigger change instead of just a patch release, as the behavior of the package actually changed. It broke the tests in some package. It is also unclear from the PR name that the behavior changed, so I suggest to rename it.)

[OlivierHnt]

I also think that the latest release should have been marked as a bigger change instead of just a patch release, as the behavior of the package actually changed. It broke the tests in some package. It is also unclear from the PR name that the behavior changed, so I suggest to rename it.

Mhm so maybe for the 0.22 releases the default should be :slow. This would be the safest option.

[lbenet]

Mhm so maybe for the 0.22 releases the default should be :slow. This would be the safest option.

I agree, the default should be :slow

[OlivierHnt]

@schillic thanks again for bringing this to our attention. The version 0.22.18 should be released soon, which introduces multiplication mode set to :slow by default. So you should get back the sharper matrix multiplication result, as before.

[schillic]

Thank you, much appreciated!