Mercurial > hg > octave-lojdl
changeset 7611:4f903c303c3c
implement subspace function
author | Jaroslav Hajek <highegg@gmail.com> |
---|---|
date | Wed, 19 Mar 2008 21:22:20 -0400 |
parents | 60398362938c |
children | c1702f963a5e |
files | scripts/ChangeLog scripts/linear-algebra/Makefile.in scripts/linear-algebra/subspace.m |
diffstat | 3 files changed, 63 insertions(+), 1 deletions(-) [+] |
line wrap: on
line diff
--- a/scripts/ChangeLog +++ b/scripts/ChangeLog @@ -1,3 +1,7 @@ +2008-03-19 Jaroslav Hajek <highegg@gmail.com> + + * linear-algebra/subspace.m: New function. + 2008-03-19 Emil Lucretiu <emil@la.mine.nu> * signal/sinetone.m: Ensure integral number of samples. @@ -35,6 +39,10 @@ * plot/__go_draw_axes__.m: Use correct symbol codes. +2008-03-17 Jaroslav Hajek <highegg@gmail.com> + + * linear-algebra/subspace.m: New function. + 2008-03-14 Kai Habel <kai.habel@gmx.de> * plot/__go_draw_axes__.m: Expicitly set gnuplot user
--- a/scripts/linear-algebra/Makefile.in +++ b/scripts/linear-algebra/Makefile.in @@ -35,7 +35,8 @@ SOURCES = __norm__.m commutation_matrix.m cond.m condest.m cross.m \ dmult.m dot.m duplication_matrix.m housh.m krylov.m krylovb.m logm.m \ - null.m onenormest.m orth.m qzhess.m rank.m rref.m trace.m vec.m vech.m + null.m onenormest.m orth.m qzhess.m rank.m rref.m subspace.m trace.m \ + vec.m vech.m DISTFILES = $(addprefix $(srcdir)/, Makefile.in $(SOURCES))
new file mode 100644 --- /dev/null +++ b/scripts/linear-algebra/subspace.m @@ -0,0 +1,53 @@ +## Copyright (C) 2008 VZLU Prague, a.s., Czech Republic +## +## This file is part of Octave. +## +## Octave is free software; you can redistribute it and/or modify it +## under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 3 of the License, or (at +## your option) any later version. +## +## Octave is distributed in the hope that it will be useful, but +## WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +## General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with Octave; see the file COPYING. If not, see +## <http://www.gnu.org/licenses/>. + +## -*- texinfo -*- +## @deftypefn {Function File} {@var{angle} =} subspace (@var{a}, @var{B}) +## Determine the largest principal angle between two subspaces +## spanned by columns of matrices @var{a} and @var{b}. +## @end deftypefn + +## Author: Jaroslav Hajek <highegg@gmail.com> + +## reference: +## [1] Andrew V. Knyazev, Merico E. Argentati: +## Principal Angles between Subspaces in an A-Based Scalar Product: +## Algorithms and Perturbation Estimates. +## SIAM Journal on Scientific Computing, Vol. 23 no. 6, pp. 2008-2040 +## +## other texts are also around... + +function ang = subspace (a, b) + + a = orth (a); + b = orth (b); + c = a'*b; + scos = min (svd (c)); + if (scos^2 > 1/2) + if (columns (a) >= columns (b)) + c = b - a*c; + else + c = a - b*c'; + endif + ssin = max (svd (c)); + ang = asin (min (ssin, 1)); + else + ang = acos (scos); + endif + +endfunction