Re: [linbox-use] Re: nullspace computation
Clement Pernet
Please answer to linbox-devel only.
> 1. Special function for nonsingular matrices.
> solveNonsingular(x, A, b) // compute x = A^{-1}b.
>
Ok why not, although I am not sure there is much efficiency to gain from
the knowledge of the non-singularity.
> 2. General case, where n is the (right) nullity of matrix/blackbox A.
>
> solve(x, A, b) // compute ARBITRARY x, such that Ax = b. If no such
> x exists, set x = 0.
>
I disagree. The value of x should not be used as an exception flag.
Let's use exceptions or dedicated flags.
I would prefer exceptions: we have them, and the cost of catching them
is really not an issue.
> Here 'arbitrary' means easiest for the implementor. This is what we
> offer now and what returns x = 0 when b = 0 and the method is
> blaselimination.
>
> solve(x, A, b, M) // user can use the method M to specify other
> handling of the inconsistent case (eg to give moore-penrose solution or
> to indicate inconsistency in the M on output - ignoring x).
>
>
Ok with that.
> nullspaceRandom(x, A) // compute a RANDOM element of nullspace(A).
> Whatever else 'random' means, it implies that, with high likelihood, n
> calls will produce a nullspace basis.
>
> nullspaceBasis(X, A) // compute a n column, full column rank
> DenseMatrix X such that AX = 0. (with BlasBlackbox becoming the
> DenseMatrix soon, and, for the case n = 0, with an X having zero columns
> being allowed by the representation!) Note that this implies a memory
> model in which the function is allowed to resize the output X. ...and
> we also give suitable remarks to the user in the comments about the case
> of a sparse matrix with large nullity ... or the case of a dense integer
> matrix with large nullity where again the basis can require n times more
> storage than A.
>
> Note that solveNonsingular allows for some extra efficiency in case
> nonsingularity is known a priori. The user's assertion of
> nonsingularity could be embedded in the method. But I think it is good
> to embed it in the function name for clarity, ease of use with default
> method, etc.
>
Non-singularity only implies simpler code ; but I doubt it enables
better efficiency. Correct me if I am wrong.
> The three general functions solve(), nullspaceVector(), nullspaceBasis()
> allow for good deployment of LinBox' various elimination based and
> blackbox methods. (all these functions take the optional method
> parameter of course, in the solutions/ directory style.) For the sparse
> case it may be good - for efficiency, including block method efficiency
> - to offer a version of nullspaceVector() that returns a clump of
> (independent) nullspace vectors.
>
> The following combination functions can then be left to the user or
> offered as well.
>
> solve(X, A, b) // compute X, a DenseMatrix of coldim n, such that Ax
> = b, for each col x of X.
> { return X = (x + y_1, ..., x + y_n),
> where the y_i are cols of nullspaceBasis(Y, A), and x is from
> solve(x, A, b) // arbitrary
> }
>
>
Yes, let's provide it.
btw, we also have to provide the solve (X,A,B) where X and B are n*k
matrices (solving for multiple columns).
> solveRandom(x, A, b) // compute a random element x of the solution
> space (an affine plane). Remark: in some cases we have a better way
> than combining an arbitrary solution and a random nullspace vector.
>
> Note that this proposal relegates some other choices
>
My proposal in the next post (to avoid spaghetti message).
Clément
David Saunders
Clement Pernet wrote:
>>
>> solve(x, A, b) // compute ARBITRARY x, such that Ax = b. If no such
>> x exists, set x = 0.
>>
> I disagree. The value of x should not be used as an exception flag.
> Let's use exceptions or dedicated flags.
> I would prefer exceptions: we have them, and the cost of catching them
> is really not an issue.
I agree that use of x to indicate no solution is bad form, but I don't
want to use exceptions in this case. Only use exceptions when the input
requirements are met, but the algorithm is not able to deliver the
output promise. And I am reluctant to say that consistency is an input
requirement to solve().
I do not want users to have to routinely use try-catch. Note that to do
so creates a burden - to write an exception catching wrapper - on
interfaces such as to SAGE. ...maybe right now there is always a
convenient wrapper [in C++] to put this in, but I hope ultimately that
many interfaces will not require code outside of their native language.
I did not address the function's output value in my previous email but
had in mind the LinBox standard: "output is a reference to the first
(i.e. output) parameter". That allows function calls to sit inside
expressions. Thus I am also reluctant to change the function value to a
boolean.
I'd be willing to choose between:
(A)
/** It is required that Az = b is consistent. If not, an exception will
be thrown. If consistency is not known a priori, either place the call
to solve() in a try-catch block or use solve(x, A, b, M) and then check
whether the Method object M indicates inconsistency */
Vector& solve(Vector& x, BB& A, const Vector& b) ...
or (B)
/// if C, then (upon return) Ax == b, else system Az = b is inconsistent.
Vector& solve(Vector& x, BB& A, const Vector& b, bool& C) ...
Note that I require that solve(x, A, b, M) [with method argument] does
NOT throw an inconsistency exception.
I'll react to Clement's function list in another email.
-dave