[37] | 1 | function [h] = householder_matrix(a, v) |
---|
| 2 | % h = householder_matrix(a, v) returns the Householder matrix that |
---|
| 3 | % will zero all elements of a except those corresponding to (any) |
---|
| 4 | % non-zero elements of v. If v has zero norm, an identity matrix is |
---|
| 5 | % returned. The matrix returned is either a left or right Householder |
---|
| 6 | % matrix dependent on whether the input parameters are column or row |
---|
| 7 | % vectors respectively. |
---|
| 8 | |
---|
| 9 | % Copyright © 2005 Stephen J. Sangwine and Nicolas Le Bihan. |
---|
| 10 | % See the file : Copyright.m for further details. |
---|
| 11 | |
---|
| 12 | % Reference: |
---|
| 13 | % |
---|
| 14 | % Sangwine, S. J. and Le Bihan, N., |
---|
| 15 | % Quaternion singular value decomposition based on bidiagonalization |
---|
| 16 | % to a real or complex matrix using quaternion Householder transformations, |
---|
| 17 | % Applied Mathematics and Computation, 182(1), 1 November 2006, 727-738, |
---|
| 18 | % DOI:10.1016/j.amc.2006.04.032. |
---|
| 19 | % |
---|
| 20 | % S. J. Sangwine and N. Le Bihan, |
---|
| 21 | % Quaternion Singular Value Decomposition based on Bidiagonalization |
---|
| 22 | % to a Real Matrix using Quaternion Householder Transformations, |
---|
| 23 | % arXiv:math.NA/0603251, 10 March 2006. Available at http://www.arxiv.org/. |
---|
| 24 | |
---|
| 25 | error(nargchk(2, 2, nargin)), error(nargoutchk(1, 1, nargout)) |
---|
| 26 | |
---|
| 27 | if size(a) ~= size(v) |
---|
| 28 | error('Input parameters must be vectors of the same size, both row, or both column') |
---|
| 29 | end |
---|
| 30 | |
---|
| 31 | [m, n] = size(a); |
---|
| 32 | |
---|
| 33 | col_vector = n == 1 & m == length(a); |
---|
| 34 | row_vector = m == 1 & n == length(a); |
---|
| 35 | |
---|
| 36 | if col_vector == 0 & row_vector == 0 |
---|
| 37 | error('Parameters must be either column or row vectors.'); |
---|
| 38 | end |
---|
| 39 | |
---|
| 40 | h = quaternion(eye(length(a))); |
---|
| 41 | |
---|
| 42 | n = norm(v); % If v has zero norm, we return an identity matrix. |
---|
| 43 | |
---|
| 44 | if n ~= 0 |
---|
| 45 | [u, zeta] = householder_vector(a, v ./ n); |
---|
| 46 | if size(u) ~= size(a) |
---|
| 47 | error('Result from function householder_vector is not of correct type (row, column)'); |
---|
| 48 | end |
---|
| 49 | if col_vector |
---|
| 50 | h = (1 ./ zeta) .* (h - u * u'); |
---|
| 51 | else |
---|
| 52 | h = (h - u.' * conj(u)) .* (1 ./ zeta); |
---|
| 53 | end |
---|
| 54 | end |
---|