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

cond(F) is not implemented for lufact return type #421

Closed
garrison opened this issue Apr 20, 2017 · 6 comments · Fixed by JuliaLang/julia#21577
Closed

cond(F) is not implemented for lufact return type #421

garrison opened this issue Apr 20, 2017 · 6 comments · Fixed by JuliaLang/julia#21577

Comments

@garrison
Copy link
Member

The documentation for lufact claims the returned factorization type implements cond(), but it does not:

julia> cond(lufact([1 0; 0 -1]))
ERROR: MethodError: no method matching cond(::Base.LinAlg.LU{Float64,Array{Float64,2}})

A similar issue exists in base/sparse/umfpack.jl.

I believe these methods should be implemented (otherwise the doc should be changed).

@garrison
Copy link
Member Author

garrison commented Apr 20, 2017

Actually, this seems to be implemented; I'm not sure why https://github.com/JuliaLang/julia/blob/8b2ee7126cc293893aae9acceb102fa766ccc6bf/base/linalg/lu.jl#L312-L314 is not being called.

EDIT: OK, the real reason is that this call does not have the default p=2, as the documentation suggests. The following works:

cond(lufact([1 0; 0 -1]), 2)

I would also expect cond(A, 2) to give the same result as cond(lufact(A), 2), but it does not:

julia> cond([1 0.2; 0 -1], 2)
1.220997512422418

julia> cond(lufact([1 0.2; 0 -1]), 2)
1.1786533995862283

Let me know if I am missing something here.

@garrison garrison changed the title cond is not implemented for lufact return type cond(F) is not implemented for lufact return type Apr 20, 2017
@andreasnoack
Copy link
Member

I think the documentation is right since cond is implemented. I'm not sure how the table in the documentation could be improved and I don't think it makes sense to implement the p=2 version for LU since such a method would have to recreate the original matrix. Right now, p=2 is wrong and returns p=Inf. Do you think the table can be improved?

@garrison
Copy link
Member Author

Right now, p=2 is wrong and returns p=Inf.

To clarify: I think p=2 is wrong (not sure what it returns, but it is not the same result as p=Inf). On the other hand, passing p=Inf does indeed give the correct result.

Do you think the table can be improved?

I think the table is fine as it is, but I think it would make sense to add a docstring for a second form of cond(), which takes an LU object as the first argument and which does not list a default p for the second argument. (The first form is here.)

@andreasnoack
Copy link
Member

On the other hand, passing p=Inf does indeed give the correct result.

This line is pretty broken. It asks LAPACK for the Inf norm condition number estimate but it passes the two norm of the input. This is why 1, 2, and Inf give different results. We should error if p != 1 && p != Inf.

but I think it would make sense to add a docstring for a second form

Sounds like the right solution.

@garrison
Copy link
Member Author

And actually, I think inv() should be added to the table as well. It is mentioned in the docstrings for bkfact, cholfact, and eigfact. But lufact, ldltfact, and qrfact currently neglect to mention it.

andreasnoack referenced this issue in JuliaLang/julia Apr 27, 2017
argument. The old version was broken and inefficient. While it would
be possible to fix just the bug, the performance problem is a consequence
of how condition number estimation works, i.e. a norm estimate is required
and that is not natural to produce from the factorized matrix. Hence,
the change to a helper function.

Fixes #21453
andreasnoack referenced this issue in JuliaLang/julia Apr 27, 2017
argument. The old version was broken and inefficient. While it would
be possible to fix just the bug, the performance problem is a consequence
of how condition number estimation works, i.e. a norm estimate is required
and that is not natural to produce from the factorized matrix. Hence,
the change to a helper function.

Fixes #21453
@andreasnoack
Copy link
Member

I decided to take a different approach. See JuliaLang/julia#21577. @garrison It would be great if you could open a PR with the doc changes.

andreasnoack referenced this issue in JuliaLang/julia Apr 27, 2017
argument. The old version was broken and inefficient. While it would
be possible to fix just the bug, the performance problem is a consequence
of how condition number estimation works, i.e. a norm estimate is required
and that is not natural to produce from the factorized matrix. Hence,
the change to a helper function.

Fixes #21453
andreasnoack referenced this issue in JuliaLang/julia May 8, 2017
argument. The old version was broken and inefficient. While it would
be possible to fix just the bug, the performance problem is a consequence
of how condition number estimation works, i.e. a norm estimate is required
and that is not natural to produce from the factorized matrix. Hence,
the change to a helper function.

Fixes #21453
andreasnoack referenced this issue in JuliaLang/julia May 8, 2017
argument. The old version was broken and inefficient. While it would
be possible to fix just the bug, the performance problem is a consequence
of how condition number estimation works, i.e. a norm estimate is required
and that is not natural to produce from the factorized matrix. Hence,
the change to a helper function.

Fixes #21453
andreasnoack referenced this issue in JuliaLang/julia May 10, 2017
…#21577)

argument. The old version was broken and inefficient. While it would
be possible to fix just the bug, the performance problem is a consequence
of how condition number estimation works, i.e. a norm estimate is required
and that is not natural to produce from the factorized matrix. Hence,
the change to a helper function.

Fixes #21453
@KristofferC KristofferC transferred this issue from JuliaLang/julia Nov 26, 2024
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 a pull request may close this issue.

2 participants