## Any subrout to compute the orthogonal complement of a space?

Open discussion regarding features, bugs, issues, vendors, etc.

### Any subrout to compute the orthogonal complement of a space?

Now, there are m (m<n) vectors in a n-dimensional cartesian space at hand, which subroutine in Lapack can give the orthogonal complement of these m vectors. Orthogonal complement mean that the dot product of any vector U in the complement set and any vector v in the m-vector set at hand is zero.

I guess there are two ways to deal with this problem:

1)solve the equation Ax=0, while A is a mxn matrix whose row vectors are the m vectors at hand. But I don't know which subrout in Lapack can solve all the solutions to these underdetermined linear equations.

2)solve the equation AX=lamda*X,A is a nxn matrix whose first m rows are the m vectors the other row vectors are 0 vectors. And then find the eigenvector corresponding to the eigenvalue lamda=0. But I think this method takes a lot of time because it computes many unnecessary eigenvalues lamda/=0.

Any more efficient way? I prefer the first method!
fatalme

Posts: 10
Joined: Sun Jul 24, 2005 10:54 am

At
http://www.netlib.org
, I found this

1. napack/null.f ( plus dependencies )
for: Compute an orthonormal basis for the space perpendicular to a given collection of vectors
Score: 100%

http://www.netlib.org/napack/null.f
fatalme

Posts: 10
Joined: Sun Jul 24, 2005 10:54 am

Another more efficient way is to perform a SVD of the mxn matrix A, then
A=USV^T, the last n-m column vectors of T forms a null space of the matrix A.
fatalme

Posts: 10
Joined: Sun Jul 24, 2005 10:54 am

Hello,

I'd say the two first methods proposed are really exotic, I do not know about napack/null.f, the last solution with the SVD seems ok. But for your problem you can more simply go with a QR factorization of the A^T.

One way to solve the problem would be to compute a QR factorization with Householder transformation of the transpose of the initial matrix. To do this you call DGEQRF on A^T where A is M-by-N where M<N.
Let B be A^T, the calling sequence should be something like:
Code: Select all
`       CALL DGEQRF( N, M, B, N, TAU, WORK, LWORK, INFO )`

Then you call DORMQR
Code: Select all
`      CALL DORMQR( 'L', 'N', N, N, M, B, N, TAU, C, N, WORK, LWORK, INFO )`

Where C has previously been initialized to the identity matrix of size N-by-N.
The first M columns of C represent an orthogonal basis for B. The next M+1 to N columns of C represent a basis of the orthogonal complement of B.

No time to fully check all the calling sequence, but that should give you a rough idea.

Julien
Julien Langou

Posts: 835
Joined: Thu Dec 09, 2004 12:32 pm
Location: Denver, CO, USA