Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

in-place minres #269

Merged
merged 20 commits into from
Apr 9, 2021
Merged

Conversation

geoffroyleconte
Copy link
Member

Copy link
Member

@amontoison amontoison left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Great PR !!!

@amontoison
Copy link
Member

if you want to rebase :

git remote add upstream https://github.com/JuliaSmoothOptimizers/Krylov.jl.git
git fetch upstream
git rebase upstream/master
git push -f

@codecov
Copy link

codecov bot commented Mar 16, 2021

Codecov Report

Merging #269 (edf0608) into master (e67df4a) will increase coverage by 0.13%.
The diff coverage is 100.00%.

Impacted file tree graph

@@            Coverage Diff             @@
##           master     #269      +/-   ##
==========================================
+ Coverage   96.96%   97.09%   +0.13%     
==========================================
  Files          32       33       +1     
  Lines        3491     3516      +25     
==========================================
+ Hits         3385     3414      +29     
+ Misses        106      102       -4     
Impacted Files Coverage Δ
src/Krylov.jl 100.00% <ø> (ø)
src/krylov_solvers.jl 100.00% <100.00%> (ø)
src/minres.jl 100.00% <100.00%> (ø)
src/diom.jl 94.25% <0.00%> (+1.14%) ⬆️
src/krylov_utils.jl 85.18% <0.00%> (+2.22%) ⬆️

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update e67df4a...edf0608. Read the comment docs.

@amontoison
Copy link
Member

@geoffroyleconte You were right, we need both S and T for KrylovSolver{S, T}. For some methods we cannot avoid the storage of vectors on CPU (DIOM, DQGMRES, ...). We can keep err_vec after all for minres but we need to check at the beginning that length(solver.err_vec) == windows.

@amontoison
Copy link
Member

@dpo It's fine for me. What do you think about it ?

- r1
- r2
- w1
- w2
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This doesn’t seem important/useful to document since the user is never supposed to touch this structure. I would remove lines 9-13 and add comments inside the struct to guide the developer.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The user may need to access x no? And maybe other parameters if we create the same struct for other Krylov solver?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

x is always returned by a Krylov method and it's already a pointer to a vector inside a KrylovSolver.

resid = norm(r) / norm(b)
@printf("MINRES: Relative residual: %8.1e\n", resid)
@test(resid ≤ minres_tol)
@test(stats.solved)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It would be great to add a test showing that minres!(A, b, solver) does not allocate.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this is not possible yet because there are still allocs done with the Residuals and Aresiduals vectors (see #272) . Even after the PR there will still be allocs because we create empty vectors for Residuals and Aresiduals and we do not fill them.

@dpo
Copy link
Member

dpo commented Mar 16, 2021

I think we need to separate cleanly the Minres storage from the stats structure. I think x should be in the stats structure so that the user never needs to look into the solver private storage. We should have an option to not accumulate the arrays of rNorms and ArNorms and be able to check that when that option is turned on, the solver doesn’t allocate. I would also put the stats structure inside the solver storage and return it at the end.

@geoffroyleconte
Copy link
Member Author

If x is in the stats structure, we will need to add it in the arguments of minres! . Do you mean something like
minres!(A, b, solver, x; ....) ?

@amontoison
Copy link
Member

Can you rebase your branch Geoffroy ?

@dpo
Copy link
Member

dpo commented Mar 16, 2021

Store the stats structure inside MinresSolver and x inside the stats. Then just return the stats instead of (x, stats).

If the operators are preallocated, it should be possible to remove all allocations.

@amontoison
Copy link
Member

amontoison commented Mar 16, 2021

@dpo It will break every code that use Krylov.jl and stats.rNorms and stats.ArNorms have a different length at each call of a Krylov method if history = true.
I don't think that it's a good idea. In the best case, we can win 160 bits...

@dpo
Copy link
Member

dpo commented Mar 16, 2021

It will break every code that use Krylov.jl

I don’t think so. If you want, you can return (stats.x, stats) for the time being. What will break? At any rate, that’s what version numbers are for. Overall I think it’s an improvement.

stats.rNorms and stats.ArNorms have a different length at each call of a Krylov method if history = true.

I think the objective of this should be to remove allocations. What I’m saying is that if history = false and the operators preallocate, we should be able to have no allocations in the solver itself.


may be used in order to create these vectors.
"""
mutable struct MinresSolver{T, S} <: KrylovSolver{T, S}
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

mutable is not required here.

@geoffroyleconte
Copy link
Member Author

@dpo , with @amontoison we think that it might be wise to keep this PR as it is for now, and to focus on removing the Residuals and Aresiduals in another PR. The ultimate goal is, as you said, to remove all allocations with minres!, maybe in version 0.8 .

@amontoison
Copy link
Member

@dpo We should think about a generic stats for all solvers. It will also solve #166.
Let's do it step by step.

src/minres.jl Outdated Show resolved Hide resolved
@dpo
Copy link
Member

dpo commented Mar 17, 2021

We should think about a generic stats for all solvers.

That’s a different topic. We’re just doing this for MINRES for now. Generic stats can come later (maybe).

@amontoison
Copy link
Member

amontoison commented Mar 17, 2021

We should think about a generic stats for all solvers.

That’s a different topic. We’re just doing this for MINRES for now. Generic stats can come later (maybe).

Yes, sure. But if you want to add x in stats like you asked, Geoffroy needs to update SimpleStats (which requires to update all the other methods that use it) or create a new stats object and implement its show function.
We should discuss about it together friday.

@dpo
Copy link
Member

dpo commented Mar 17, 2021

Sure, let’s duplicate SimpleStats to make one just for MINRES for now. This is our proof of concept. It’s not that difficult. No need to implement show at this point. The main point is: can we get rid of all allocations in the second call?

@amontoison
Copy link
Member

@amontoison
Copy link
Member

amontoison commented Mar 17, 2021

Sure, let’s duplicate SimpleStats to make one just for MINRES for now. This is our proof of concept. It’s not that difficult. No need to implement show at this point. The main point is: can we get rid of all allocations in the second call?

Yes, it will be the case for MINRES (and all methods that do not require A'). We just need to pass a PreallocatedLinearOperator otherwise a new one will be created at each call of minres!.

We should create the PrellocatedLinearOperator with MinresSolver(A, b) to avoid that.

@dpo
Copy link
Member

dpo commented Mar 17, 2021

Let’s do it so @geoffroyleconte can use it in his code.

@amontoison
Copy link
Member

@geoffroyleconte
You can continue without me.

@github-actions
Copy link
Contributor

github-actions bot commented Apr 7, 2021

Package name latest stable
CaNNOLeS.jl
DCI.jl
JSOSolvers.jl
Percival.jl

@github-actions
Copy link
Contributor

github-actions bot commented Apr 7, 2021

Package name latest stable
CaNNOLeS.jl
DCI.jl
JSOSolvers.jl
Percival.jl

@geoffroyleconte
Copy link
Member Author

@dpo @amontoison can we merge this if I check that the Julia version is 1.5 for the allocs test with minres! ? It is very inconvenient to have to rebase every time... We can open an issue to make allocs tests for other versions later.

@dpo
Copy link
Member

dpo commented Apr 7, 2021

I don't understand. What are the differences with what I proposed and I did for BICGSTAB ?

I didn’t realize that you simply deactivated the test for Julia ≤ 1.5. That’s ok I guess.

@dpo
Copy link
Member

dpo commented Apr 7, 2021

Also it looks like the coverage went down. I think we just need to call the variants in the tests: A = dense matrix, b = sparse vector, etc.

@dpo
Copy link
Member

dpo commented Apr 7, 2021

@geoffroyleconte Could you deactivate the MINRES allocation test for Julia ≤ 1.5 as @amontoison did? The tests should pass with Julia 1.3.

@github-actions
Copy link
Contributor

github-actions bot commented Apr 8, 2021

Package name latest stable
CaNNOLeS.jl
DCI.jl
JSOSolvers.jl
Percival.jl

solver = eval(Symbol(uppercasefirst(string(fn))[1:end-1], "Solver"))(A, b)
@eval $fn($solver, $A, $b)
@eval $fn($solver, $transpose($A), $b)
@eval $fn($solver, $adjoint($A), $b)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

See the annotations in variants.jl to see which are not covered by the current tests.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There is an issue when b is a SparseVector, I did not find how to fix it. But anyway the function converts bto a Vector so this is not really an issue, the user just has to give a Vectorinstead of a SparseVector. All these issues will be fixed when we change LinearOperators

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I understand, but what we’re trying to fix now is coverage. What’s the problem with SparseVector?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Since we have to give a PreallocatedLinearOperator that has a preallocated vector that we can keep track of so that the in-place version works correctly, I suggest we remove the variants.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Coverage is completely separate from the in-place solver. We need to test the variants simply to ensure that they are callable and behave correctly. It's understood that they will allocate. But allocations are not what we're testing. If you're having an issue with SparseVectors, then it's an indication that there's a bug and we should fix it.

Could you paste the error here?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There is no issue with minres, only with minres!. It looks like the coverage is increased. We are loosing time on minor issues that will be fixed with the changes that will be made on LinearOperators.jl. We should say that for now the user has to use PreallocatedLinearOperators. It is actually the only way he can avoid allocations with the current implementation. I can add this to the documentation if needed.

Here is the error:

minres! ERROR: LoadError: LoadError: TypeError: in typeassert, expected SparseVector{Float32,Int32}, got a value of type Array{Float32,1} Stacktrace: [1] *(::PreallocatedLinearOperator{Float32}, ::SparseVector{Float32,Int32}) at C:\Users\Geoffroy Leconte\.julia\packages\LinearOperators\qFko7\src\operations.jl:7 [2] minres!(::MinresSolver{Float32,SparseVector{Float32,Int32}}, ::PreallocatedLinearOperator{Float32}, ::Array{Float32,1}; M::opEye, λ::Float32, atol::Float32, rtol::Float32, ratol

@github-actions
Copy link
Contributor

github-actions bot commented Apr 8, 2021

Package name latest stable
CaNNOLeS.jl
DCI.jl
JSOSolvers.jl
Percival.jl

@github-actions
Copy link
Contributor

github-actions bot commented Apr 9, 2021

Package name latest stable
CaNNOLeS.jl
DCI.jl
JSOSolvers.jl
Percival.jl

@geoffroyleconte
Copy link
Member Author

I added annotations so that the user is forced to use PreallocatedLinearOperators with minres!.

I can raise an issue to deal with b when it is a SparseVector. This may be solved with JuliaSmoothOptimizers/LinearOperators.jl#172 , but if it works it will change all the Krylov variants with b::SparseVector{T} so we should deal with this in another PR.

@@ -22,7 +23,7 @@ mutable struct MinresSolver{T, S} <: KrylovSolver{T, S}
stats :: SimpleStats{T}
end

function MinresSolver(A, b :: AbstractVector{T}; window :: Int=5) where T <: AbstractFloat
function MinresSolver(A :: PreallocatedLinearOperator{T}, b :: Vector{T}; window :: Int=5) where T <: AbstractFloat
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We cannot use anymore minres on GPU :(

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I thought the annotation should go in minres!(), not here. Would that also break the GPU?

Copy link
Member

@amontoison amontoison Apr 9, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

b must be an AbstractVector{T} to keep support with CuVector{T}. (GPU storage).
The parametrization of A is not needed too : https://juliasmoothoptimizers.github.io/Krylov.jl/dev/matrix-free/

But one thing should be explained somewhere, if the user want to use minres!, he should wrap itselft A::AbstractMatrix in a PreallocatedLinearOperator. It's only done automatically with minres.
If the user already uses its linear operator (LinearOperators.jl, LinearMaps.jl, ...) he doesn't need to do anything.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OK, I removed annotations and documented that the user has to use PreallocatedLinearOperators to avoid allocations.

@github-actions
Copy link
Contributor

github-actions bot commented Apr 9, 2021

Package name latest stable
CaNNOLeS.jl
DCI.jl
JSOSolvers.jl
Percival.jl

@dpo dpo merged commit 0c077e0 into JuliaSmoothOptimizers:master Apr 9, 2021
@dpo
Copy link
Member

dpo commented Apr 9, 2021

Great! Thank you!

@geoffroyleconte geoffroyleconte deleted the minres-ip branch April 9, 2021 19:35
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants