git.maemo.org Git - opencv/blob - 3rdparty/lapack/sorm2l.c

   1 #include "clapack.h"
   2
   3 /* Table of constant values */
   4
   5 static integer c__1 = 1;
   6
   7 /* Subroutine */ int sorm2l_(char *side, char *trans, integer *m, integer *n,
   8         integer *k, real *a, integer *lda, real *tau, real *c__, integer *ldc,
   9          real *work, integer *info)
  10 {
  11     /* System generated locals */
  12     integer a_dim1, a_offset, c_dim1, c_offset, i__1, i__2;
  13
  14     /* Local variables */
  15     integer i__, i1, i2, i3, mi, ni, nq;
  16     real aii;
  17     logical left;
  18     extern logical lsame_(char *, char *);
  19     extern /* Subroutine */ int slarf_(char *, integer *, integer *, real *,
  20             integer *, real *, real *, integer *, real *), xerbla_(
  21             char *, integer *);
  22     logical notran;
  23
  24
  25 /*  -- LAPACK routine (version 3.1) -- */
  26 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
  27 /*     November 2006 */
  28
  29 /*     .. Scalar Arguments .. */
  30 /*     .. */
  31 /*     .. Array Arguments .. */
  32 /*     .. */
  33
  34 /*  Purpose */
  35 /*  ======= */
  36
  37 /*  SORM2L overwrites the general real m by n matrix C with */
  38
  39 /*        Q * C  if SIDE = 'L' and TRANS = 'N', or */
  40
  41 /*        Q'* C  if SIDE = 'L' and TRANS = 'T', or */
  42
  43 /*        C * Q  if SIDE = 'R' and TRANS = 'N', or */
  44
  45 /*        C * Q' if SIDE = 'R' and TRANS = 'T', */
  46
  47 /*  where Q is a real orthogonal matrix defined as the product of k */
  48 /*  elementary reflectors */
  49
  50 /*        Q = H(k) . . . H(2) H(1) */
  51
  52 /*  as returned by SGEQLF. Q is of order m if SIDE = 'L' and of order n */
  53 /*  if SIDE = 'R'. */
  54
  55 /*  Arguments */
  56 /*  ========= */
  57
  58 /*  SIDE    (input) CHARACTER*1 */
  59 /*          = 'L': apply Q or Q' from the Left */
  60 /*          = 'R': apply Q or Q' from the Right */
  61
  62 /*  TRANS   (input) CHARACTER*1 */
  63 /*          = 'N': apply Q  (No transpose) */
  64 /*          = 'T': apply Q' (Transpose) */
  65
  66 /*  M       (input) INTEGER */
  67 /*          The number of rows of the matrix C. M >= 0. */
  68
  69 /*  N       (input) INTEGER */
  70 /*          The number of columns of the matrix C. N >= 0. */
  71
  72 /*  K       (input) INTEGER */
  73 /*          The number of elementary reflectors whose product defines */
  74 /*          the matrix Q. */
  75 /*          If SIDE = 'L', M >= K >= 0; */
  76 /*          if SIDE = 'R', N >= K >= 0. */
  77
  78 /*  A       (input) REAL array, dimension (LDA,K) */
  79 /*          The i-th column must contain the vector which defines the */
  80 /*          elementary reflector H(i), for i = 1,2,...,k, as returned by */
  81 /*          SGEQLF in the last k columns of its array argument A. */
  82 /*          A is modified by the routine but restored on exit. */
  83
  84 /*  LDA     (input) INTEGER */
  85 /*          The leading dimension of the array A. */
  86 /*          If SIDE = 'L', LDA >= max(1,M); */
  87 /*          if SIDE = 'R', LDA >= max(1,N). */
  88
  89 /*  TAU     (input) REAL array, dimension (K) */
  90 /*          TAU(i) must contain the scalar factor of the elementary */
  91 /*          reflector H(i), as returned by SGEQLF. */
  92
  93 /*  C       (input/output) REAL array, dimension (LDC,N) */
  94 /*          On entry, the m by n matrix C. */
  95 /*          On exit, C is overwritten by Q*C or Q'*C or C*Q' or C*Q. */
  96
  97 /*  LDC     (input) INTEGER */
  98 /*          The leading dimension of the array C. LDC >= max(1,M). */
  99
 100 /*  WORK    (workspace) REAL array, dimension */
 101 /*                                   (N) if SIDE = 'L', */
 102 /*                                   (M) if SIDE = 'R' */
 103
 104 /*  INFO    (output) INTEGER */
 105 /*          = 0: successful exit */
 106 /*          < 0: if INFO = -i, the i-th argument had an illegal value */
 107
 108 /*  ===================================================================== */
 109
 110 /*     .. Parameters .. */
 111 /*     .. */
 112 /*     .. Local Scalars .. */
 113 /*     .. */
 114 /*     .. External Functions .. */
 115 /*     .. */
 116 /*     .. External Subroutines .. */
 117 /*     .. */
 118 /*     .. Intrinsic Functions .. */
 119 /*     .. */
 120 /*     .. Executable Statements .. */
 121
 122 /*     Test the input arguments */
 123
 124     /* Parameter adjustments */
 125     a_dim1 = *lda;
 126     a_offset = 1 + a_dim1;
 127     a -= a_offset;
 128     --tau;
 129     c_dim1 = *ldc;
 130     c_offset = 1 + c_dim1;
 131     c__ -= c_offset;
 132     --work;
 133
 134     /* Function Body */
 135     *info = 0;
 136     left = lsame_(side, "L");
 137     notran = lsame_(trans, "N");
 138
 139 /*     NQ is the order of Q */
 140
 141     if (left) {
 142         nq = *m;
 143     } else {
 144         nq = *n;
 145     }
 146     if (! left && ! lsame_(side, "R")) {
 147         *info = -1;
 148     } else if (! notran && ! lsame_(trans, "T")) {
 149         *info = -2;
 150     } else if (*m < 0) {
 151         *info = -3;
 152     } else if (*n < 0) {
 153         *info = -4;
 154     } else if (*k < 0 || *k > nq) {
 155         *info = -5;
 156     } else if (*lda < max(1,nq)) {
 157         *info = -7;
 158     } else if (*ldc < max(1,*m)) {
 159         *info = -10;
 160     }
 161     if (*info != 0) {
 162         i__1 = -(*info);
 163         xerbla_("SORM2L", &i__1);
 164         return 0;
 165     }
 166
 167 /*     Quick return if possible */
 168
 169     if (*m == 0 || *n == 0 || *k == 0) {
 170         return 0;
 171     }
 172
 173     if (left && notran || ! left && ! notran) {
 174         i1 = 1;
 175         i2 = *k;
 176         i3 = 1;
 177     } else {
 178         i1 = *k;
 179         i2 = 1;
 180         i3 = -1;
 181     }
 182
 183     if (left) {
 184         ni = *n;
 185     } else {
 186         mi = *m;
 187     }
 188
 189     i__1 = i2;
 190     i__2 = i3;
 191     for (i__ = i1; i__2 < 0 ? i__ >= i__1 : i__ <= i__1; i__ += i__2) {
 192         if (left) {
 193
 194 /*           H(i) is applied to C(1:m-k+i,1:n) */
 195
 196             mi = *m - *k + i__;
 197         } else {
 198
 199 /*           H(i) is applied to C(1:m,1:n-k+i) */
 200
 201             ni = *n - *k + i__;
 202         }
 203
 204 /*        Apply H(i) */
 205
 206         aii = a[nq - *k + i__ + i__ * a_dim1];
 207         a[nq - *k + i__ + i__ * a_dim1] = 1.f;
 208         slarf_(side, &mi, &ni, &a[i__ * a_dim1 + 1], &c__1, &tau[i__], &c__[
 209                 c_offset], ldc, &work[1]);
 210         a[nq - *k + i__ + i__ * a_dim1] = aii;
 211 /* L10: */
 212     }
 213     return 0;
 214
 215 /*     End of SORM2L */
 216
 217 } /* sorm2l_ */