Mercurial > octave-nkf
changeset 740:d8295febb0df
[project @ 1994-09-30 14:42:37 by jwe]
| author | jwe |
|---|---|
| date | Fri, 30 Sep 1994 14:42:37 +0000 |
| parents | d452e6f52d12 |
| children | 2d2f3c07cedd |
| files | liboctave/CMatrix.cc liboctave/CMatrix.h liboctave/dMatrix.cc liboctave/dMatrix.h |
| diffstat | 4 files changed, 82 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
--- a/liboctave/CMatrix.cc Mon Sep 26 17:59:45 1994 +0000 +++ b/liboctave/CMatrix.cc Fri Sep 30 14:42:37 1994 +0000 @@ -31,11 +31,13 @@ #include <sys/types.h> #include <iostream.h> +#include <float.h> #include <Complex.h> #include "mx-base.h" #include "CmplxDET.h" +#include "CmplxSVD.h" #include "mx-inlines.cc" #include "lo-error.h" #include "f77-uscore.h" @@ -957,6 +959,43 @@ } ComplexMatrix +ComplexMatrix::pseudo_inverse (double tol) +{ + ComplexSVD result (*this); + + DiagMatrix S = result.singular_values (); + ComplexMatrix U = result.left_singular_matrix (); + ComplexMatrix V = result.right_singular_matrix (); + + ColumnVector sigma = S.diag (); + + int r = sigma.length () - 1; + int nr = rows (); + int nc = cols (); + + if (tol <= 0.0) + { + if (nr > nc) + tol = nr * sigma.elem (0) * DBL_EPSILON; + else + tol = nc * sigma.elem (0) * DBL_EPSILON; + } + + while (r >= 0 && sigma.elem (r) < tol) + r--; + + if (r < 0) + return ComplexMatrix (nc, nr, 0.0); + else + { + ComplexMatrix Ur = U.extract (0, 0, nr-1, r); + DiagMatrix D = DiagMatrix (sigma.extract (0, r)) . inverse (); + ComplexMatrix Vr = V.extract (0, 0, nc-1, r); + return Vr * D * Ur.hermitian (); + } +} + +ComplexMatrix ComplexMatrix::fourier (void) const { int nr = rows ();
--- a/liboctave/CMatrix.h Mon Sep 26 17:59:45 1994 +0000 +++ b/liboctave/CMatrix.h Fri Sep 30 14:42:37 1994 +0000 @@ -134,6 +134,8 @@ ComplexMatrix inverse (int& info) const; ComplexMatrix inverse (int& info, double& rcond) const; + ComplexMatrix pseudo_inverse (double tol = 0.0); + ComplexMatrix fourier (void) const; ComplexMatrix ifourier (void) const;
--- a/liboctave/dMatrix.cc Mon Sep 26 17:59:45 1994 +0000 +++ b/liboctave/dMatrix.cc Fri Sep 30 14:42:37 1994 +0000 @@ -32,11 +32,13 @@ #include <sys/types.h> #include <iostream.h> #include <stdio.h> +#include <float.h> #include <Complex.h> #include "mx-base.h" #include "dbleDET.h" +#include "dbleSVD.h" #include "mx-inlines.cc" #include "lo-error.h" #include "f77-uscore.h" @@ -605,6 +607,43 @@ return Matrix (tmp_data, nr, nc); } +Matrix +Matrix::pseudo_inverse (double tol) +{ + SVD result (*this); + + DiagMatrix S = result.singular_values (); + Matrix U = result.left_singular_matrix (); + Matrix V = result.right_singular_matrix (); + + ColumnVector sigma = S.diag (); + + int r = sigma.length () - 1; + int nr = rows (); + int nc = cols (); + + if (tol <= 0.0) + { + if (nr > nc) + tol = nr * sigma.elem (0) * DBL_EPSILON; + else + tol = nc * sigma.elem (0) * DBL_EPSILON; + } + + while (r >= 0 && sigma.elem (r) < tol) + r--; + + if (r < 0) + return Matrix (nc, nr, 0.0); + else + { + Matrix Ur = U.extract (0, 0, nr-1, r); + DiagMatrix D = DiagMatrix (sigma.extract (0, r)) . inverse (); + Matrix Vr = V.extract (0, 0, nc-1, r); + return Vr * D * Ur.transpose (); + } +} + ComplexMatrix Matrix::fourier (void) const {
--- a/liboctave/dMatrix.h Mon Sep 26 17:59:45 1994 +0000 +++ b/liboctave/dMatrix.h Fri Sep 30 14:42:37 1994 +0000 @@ -109,6 +109,8 @@ Matrix inverse (int& info) const; Matrix inverse (int& info, double& rcond) const; + Matrix pseudo_inverse (double tol = 0.0); + ComplexMatrix fourier (void) const; ComplexMatrix ifourier (void) const;