qr.h 3.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475
  1. /*Daala video codec
  2. Copyright (c) 2007-2014 Daala project contributors. All rights reserved.
  3. Redistribution and use in source and binary forms, with or without
  4. modification, are permitted provided that the following conditions are met:
  5. - Redistributions of source code must retain the above copyright notice, this
  6. list of conditions and the following disclaimer.
  7. - Redistributions in binary form must reproduce the above copyright notice,
  8. this list of conditions and the following disclaimer in the documentation
  9. and/or other materials provided with the distribution.
  10. THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS “AS IS”
  11. AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  12. IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
  13. DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE
  14. FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
  15. DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
  16. SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
  17. CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
  18. OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
  19. OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.*/
  20. #if !defined(_qr_H)
  21. # define _qr_H (1)
  22. # include "matidx.h"
  23. /*Computes the (skinny) QR-decomposition A=Q.R of an m by n matrix A using
  24. Householder reflections.
  25. The m by MIN(m, n) matrix Q is column-orthonormal (so that Q^T.Q=I), and
  26. the MIN(m, n) by n matrix R is upper-trapezoidal (so that R_{i, j} = 0 for
  27. all j < i).
  28. Because Householder reflections operate on columns, the decomposition is done
  29. on the transpose of the matrix (e.g., column-major order) for cache
  30. efficiency.
  31. The result of the computation is stored in place, and not in direct form, but
  32. individual components can be extracted with qrdecomp_hh_get_r() and
  33. qrdecomp_hh_get_q().
  34. aat: The n by m transpose of the matrix to decompose.
  35. The entries to the left of the diagonal are replaced with the entries of
  36. R^T, and the entries to the right of the diagonal are replaced by the
  37. Householder vectors used to compute Q.
  38. aat_stride: The row stride of aat (column stride of A).
  39. d: Scratch space for the diagonal of R.
  40. This must have storage for MIN(m, n) entries.
  41. qqt: Outputs the MIN(n, m) by m orthogonal matrix Q^T.
  42. This may be NULL to skip computation of this matrix if it is not needed.
  43. qqt_stride: The row stride of qqt (column stride of Q).
  44. rr: Outputs the upper-trapezoidal matrix R.
  45. This may be NULL to skip the output of this matrix in upper-triangular
  46. form if it is not needed, otherwise it must have storage for
  47. UT_SZ(MIN(m, n), n) entries and rr[UT_IDX(i, j, n)] will contain
  48. R_{i, j} for i <= j on exit.
  49. Return: The rank of the matrix aat, computed as the number of non-zero
  50. entries on the diagonal of R.*/
  51. int qrdecomp_hh(double *aat,int aat_stride, double *d,
  52. double *qqt, int qqt_stride, double *rr, int n, int m);
  53. /*Finds the least-squares solution of an over-determined linear system
  54. AX ~= B using a previously computed full-rank QR-decomposition of A.
  55. The solution for l different column vectors of B is found simultaneously.
  56. qrt: The QR-decomposition, as computed by qrdecomp_hh().
  57. The rank of the decomposition (as returned by qrdecomp_hh()) must be n.
  58. qrt_stride: The row stride of qrt.
  59. d: The diagonal of R, as computed by qrdecomp_hh().
  60. bbt: The l by m transpose of the l right hand sides B.
  61. These are stored as row vectors to improve cache efficiency and to make
  62. the common case of a single solution vector simpler to handle.
  63. They are replaced by the solution vectors on output.*/
  64. void qrsolve_hh(double *qrt, int qrt_stride, double *d,
  65. double *bbt, int bbt_stride, int n, int m, int l);
  66. #endif