dormtr.c 8.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296
  1. /* dormtr.f -- translated by f2c (version 20061008).
  2. You must link the resulting object file with libf2c:
  3. on Microsoft Windows system, link with libf2c.lib;
  4. on Linux or Unix systems, link with .../path/to/libf2c.a -lm
  5. or, if you install libf2c.a in a standard place, with -lf2c -lm
  6. -- in that order, at the end of the command line, as in
  7. cc *.o -lf2c -lm
  8. Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
  9. http://www.netlib.org/f2c/libf2c.zip
  10. */
  11. #include "f2c.h"
  12. #include "blaswrap.h"
  13. /* Table of constant values */
  14. static integer c__1 = 1;
  15. static integer c_n1 = -1;
  16. static integer c__2 = 2;
  17. /* Subroutine */ int _starpu_dormtr_(char *side, char *uplo, char *trans, integer *m,
  18. integer *n, doublereal *a, integer *lda, doublereal *tau, doublereal *
  19. c__, integer *ldc, doublereal *work, integer *lwork, integer *info)
  20. {
  21. /* System generated locals */
  22. address a__1[2];
  23. integer a_dim1, a_offset, c_dim1, c_offset, i__1[2], i__2, i__3;
  24. char ch__1[2];
  25. /* Builtin functions */
  26. /* Subroutine */ int s_cat(char *, char **, integer *, integer *, ftnlen);
  27. /* Local variables */
  28. integer i1, i2, nb, mi, ni, nq, nw;
  29. logical left;
  30. extern logical _starpu_lsame_(char *, char *);
  31. integer iinfo;
  32. logical upper;
  33. extern /* Subroutine */ int _starpu_xerbla_(char *, integer *);
  34. extern integer _starpu_ilaenv_(integer *, char *, char *, integer *, integer *,
  35. integer *, integer *);
  36. extern /* Subroutine */ int _starpu_dormql_(char *, char *, integer *, integer *,
  37. integer *, doublereal *, integer *, doublereal *, doublereal *,
  38. integer *, doublereal *, integer *, integer *),
  39. _starpu_dormqr_(char *, char *, integer *, integer *, integer *,
  40. doublereal *, integer *, doublereal *, doublereal *, integer *,
  41. doublereal *, integer *, integer *);
  42. integer lwkopt;
  43. logical lquery;
  44. /* -- LAPACK routine (version 3.2) -- */
  45. /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
  46. /* November 2006 */
  47. /* .. Scalar Arguments .. */
  48. /* .. */
  49. /* .. Array Arguments .. */
  50. /* .. */
  51. /* Purpose */
  52. /* ======= */
  53. /* DORMTR overwrites the general real M-by-N matrix C with */
  54. /* SIDE = 'L' SIDE = 'R' */
  55. /* TRANS = 'N': Q * C C * Q */
  56. /* TRANS = 'T': Q**T * C C * Q**T */
  57. /* where Q is a real orthogonal matrix of order nq, with nq = m if */
  58. /* SIDE = 'L' and nq = n if SIDE = 'R'. Q is defined as the product of */
  59. /* nq-1 elementary reflectors, as returned by DSYTRD: */
  60. /* if UPLO = 'U', Q = H(nq-1) . . . H(2) H(1); */
  61. /* if UPLO = 'L', Q = H(1) H(2) . . . H(nq-1). */
  62. /* Arguments */
  63. /* ========= */
  64. /* SIDE (input) CHARACTER*1 */
  65. /* = 'L': apply Q or Q**T from the Left; */
  66. /* = 'R': apply Q or Q**T from the Right. */
  67. /* UPLO (input) CHARACTER*1 */
  68. /* = 'U': Upper triangle of A contains elementary reflectors */
  69. /* from DSYTRD; */
  70. /* = 'L': Lower triangle of A contains elementary reflectors */
  71. /* from DSYTRD. */
  72. /* TRANS (input) CHARACTER*1 */
  73. /* = 'N': No transpose, apply Q; */
  74. /* = 'T': Transpose, apply Q**T. */
  75. /* M (input) INTEGER */
  76. /* The number of rows of the matrix C. M >= 0. */
  77. /* N (input) INTEGER */
  78. /* The number of columns of the matrix C. N >= 0. */
  79. /* A (input) DOUBLE PRECISION array, dimension */
  80. /* (LDA,M) if SIDE = 'L' */
  81. /* (LDA,N) if SIDE = 'R' */
  82. /* The vectors which define the elementary reflectors, as */
  83. /* returned by DSYTRD. */
  84. /* LDA (input) INTEGER */
  85. /* The leading dimension of the array A. */
  86. /* LDA >= max(1,M) if SIDE = 'L'; LDA >= max(1,N) if SIDE = 'R'. */
  87. /* TAU (input) DOUBLE PRECISION array, dimension */
  88. /* (M-1) if SIDE = 'L' */
  89. /* (N-1) if SIDE = 'R' */
  90. /* TAU(i) must contain the scalar factor of the elementary */
  91. /* reflector H(i), as returned by DSYTRD. */
  92. /* C (input/output) DOUBLE PRECISION array, dimension (LDC,N) */
  93. /* On entry, the M-by-N matrix C. */
  94. /* On exit, C is overwritten by Q*C or Q**T*C or C*Q**T or C*Q. */
  95. /* LDC (input) INTEGER */
  96. /* The leading dimension of the array C. LDC >= max(1,M). */
  97. /* WORK (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LWORK)) */
  98. /* On exit, if INFO = 0, WORK(1) returns the optimal LWORK. */
  99. /* LWORK (input) INTEGER */
  100. /* The dimension of the array WORK. */
  101. /* If SIDE = 'L', LWORK >= max(1,N); */
  102. /* if SIDE = 'R', LWORK >= max(1,M). */
  103. /* For optimum performance LWORK >= N*NB if SIDE = 'L', and */
  104. /* LWORK >= M*NB if SIDE = 'R', where NB is the optimal */
  105. /* blocksize. */
  106. /* If LWORK = -1, then a workspace query is assumed; the routine */
  107. /* only calculates the optimal size of the WORK array, returns */
  108. /* this value as the first entry of the WORK array, and no error */
  109. /* message related to LWORK is issued by XERBLA. */
  110. /* INFO (output) INTEGER */
  111. /* = 0: successful exit */
  112. /* < 0: if INFO = -i, the i-th argument had an illegal value */
  113. /* ===================================================================== */
  114. /* .. Local Scalars .. */
  115. /* .. */
  116. /* .. External Functions .. */
  117. /* .. */
  118. /* .. External Subroutines .. */
  119. /* .. */
  120. /* .. Intrinsic Functions .. */
  121. /* .. */
  122. /* .. Executable Statements .. */
  123. /* Test the input arguments */
  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. /* Function Body */
  134. *info = 0;
  135. left = _starpu_lsame_(side, "L");
  136. upper = _starpu_lsame_(uplo, "U");
  137. lquery = *lwork == -1;
  138. /* NQ is the order of Q and NW is the minimum dimension of WORK */
  139. if (left) {
  140. nq = *m;
  141. nw = *n;
  142. } else {
  143. nq = *n;
  144. nw = *m;
  145. }
  146. if (! left && ! _starpu_lsame_(side, "R")) {
  147. *info = -1;
  148. } else if (! upper && ! _starpu_lsame_(uplo, "L")) {
  149. *info = -2;
  150. } else if (! _starpu_lsame_(trans, "N") && ! _starpu_lsame_(trans,
  151. "T")) {
  152. *info = -3;
  153. } else if (*m < 0) {
  154. *info = -4;
  155. } else if (*n < 0) {
  156. *info = -5;
  157. } else if (*lda < max(1,nq)) {
  158. *info = -7;
  159. } else if (*ldc < max(1,*m)) {
  160. *info = -10;
  161. } else if (*lwork < max(1,nw) && ! lquery) {
  162. *info = -12;
  163. }
  164. if (*info == 0) {
  165. if (upper) {
  166. if (left) {
  167. /* Writing concatenation */
  168. i__1[0] = 1, a__1[0] = side;
  169. i__1[1] = 1, a__1[1] = trans;
  170. s_cat(ch__1, a__1, i__1, &c__2, (ftnlen)2);
  171. i__2 = *m - 1;
  172. i__3 = *m - 1;
  173. nb = _starpu_ilaenv_(&c__1, "DORMQL", ch__1, &i__2, n, &i__3, &c_n1);
  174. } else {
  175. /* Writing concatenation */
  176. i__1[0] = 1, a__1[0] = side;
  177. i__1[1] = 1, a__1[1] = trans;
  178. s_cat(ch__1, a__1, i__1, &c__2, (ftnlen)2);
  179. i__2 = *n - 1;
  180. i__3 = *n - 1;
  181. nb = _starpu_ilaenv_(&c__1, "DORMQL", ch__1, m, &i__2, &i__3, &c_n1);
  182. }
  183. } else {
  184. if (left) {
  185. /* Writing concatenation */
  186. i__1[0] = 1, a__1[0] = side;
  187. i__1[1] = 1, a__1[1] = trans;
  188. s_cat(ch__1, a__1, i__1, &c__2, (ftnlen)2);
  189. i__2 = *m - 1;
  190. i__3 = *m - 1;
  191. nb = _starpu_ilaenv_(&c__1, "DORMQR", ch__1, &i__2, n, &i__3, &c_n1);
  192. } else {
  193. /* Writing concatenation */
  194. i__1[0] = 1, a__1[0] = side;
  195. i__1[1] = 1, a__1[1] = trans;
  196. s_cat(ch__1, a__1, i__1, &c__2, (ftnlen)2);
  197. i__2 = *n - 1;
  198. i__3 = *n - 1;
  199. nb = _starpu_ilaenv_(&c__1, "DORMQR", ch__1, m, &i__2, &i__3, &c_n1);
  200. }
  201. }
  202. lwkopt = max(1,nw) * nb;
  203. work[1] = (doublereal) lwkopt;
  204. }
  205. if (*info != 0) {
  206. i__2 = -(*info);
  207. _starpu_xerbla_("DORMTR", &i__2);
  208. return 0;
  209. } else if (lquery) {
  210. return 0;
  211. }
  212. /* Quick return if possible */
  213. if (*m == 0 || *n == 0 || nq == 1) {
  214. work[1] = 1.;
  215. return 0;
  216. }
  217. if (left) {
  218. mi = *m - 1;
  219. ni = *n;
  220. } else {
  221. mi = *m;
  222. ni = *n - 1;
  223. }
  224. if (upper) {
  225. /* Q was determined by a call to DSYTRD with UPLO = 'U' */
  226. i__2 = nq - 1;
  227. _starpu_dormql_(side, trans, &mi, &ni, &i__2, &a[(a_dim1 << 1) + 1], lda, &
  228. tau[1], &c__[c_offset], ldc, &work[1], lwork, &iinfo);
  229. } else {
  230. /* Q was determined by a call to DSYTRD with UPLO = 'L' */
  231. if (left) {
  232. i1 = 2;
  233. i2 = 1;
  234. } else {
  235. i1 = 1;
  236. i2 = 2;
  237. }
  238. i__2 = nq - 1;
  239. _starpu_dormqr_(side, trans, &mi, &ni, &i__2, &a[a_dim1 + 2], lda, &tau[1], &
  240. c__[i1 + i2 * c_dim1], ldc, &work[1], lwork, &iinfo);
  241. }
  242. work[1] = (doublereal) lwkopt;
  243. return 0;
  244. /* End of DORMTR */
  245. } /* _starpu_dormtr_ */