diff options
author | shmel1k <shmel1k@ydb.tech> | 2022-09-02 12:44:59 +0300 |
---|---|---|
committer | shmel1k <shmel1k@ydb.tech> | 2022-09-02 12:44:59 +0300 |
commit | 90d450f74722da7859d6f510a869f6c6908fd12f (patch) | |
tree | 538c718dedc76cdfe37ad6d01ff250dd930d9278 /contrib/libs/clapack/dormlq.c | |
parent | 01f64c1ecd0d4ffa9e3a74478335f1745f26cc75 (diff) | |
download | ydb-90d450f74722da7859d6f510a869f6c6908fd12f.tar.gz |
[] add metering mode to CLI
Diffstat (limited to 'contrib/libs/clapack/dormlq.c')
-rw-r--r-- | contrib/libs/clapack/dormlq.c | 335 |
1 files changed, 335 insertions, 0 deletions
diff --git a/contrib/libs/clapack/dormlq.c b/contrib/libs/clapack/dormlq.c new file mode 100644 index 0000000000..bcc1cc4700 --- /dev/null +++ b/contrib/libs/clapack/dormlq.c @@ -0,0 +1,335 @@ +/* dormlq.f -- translated by f2c (version 20061008). + You must link the resulting object file with libf2c: + on Microsoft Windows system, link with libf2c.lib; + on Linux or Unix systems, link with .../path/to/libf2c.a -lm + or, if you install libf2c.a in a standard place, with -lf2c -lm + -- in that order, at the end of the command line, as in + cc *.o -lf2c -lm + Source for libf2c is in /netlib/f2c/libf2c.zip, e.g., + + http://www.netlib.org/f2c/libf2c.zip +*/ + +#include "f2c.h" +#include "blaswrap.h" + +/* Table of constant values */ + +static integer c__1 = 1; +static integer c_n1 = -1; +static integer c__2 = 2; +static integer c__65 = 65; + +/* Subroutine */ int dormlq_(char *side, char *trans, integer *m, integer *n, + integer *k, doublereal *a, integer *lda, doublereal *tau, doublereal * + c__, integer *ldc, doublereal *work, integer *lwork, integer *info) +{ + /* System generated locals */ + address a__1[2]; + integer a_dim1, a_offset, c_dim1, c_offset, i__1, i__2, i__3[2], i__4, + i__5; + char ch__1[3]; + ch__1[2] = 0; + + /* Builtin functions */ + /* Subroutine */ int s_cat(char *, char **, integer *, integer *, ftnlen); + + /* Local variables */ + integer i__; + doublereal t[4160] /* was [65][64] */; + integer i1, i2, i3, ib, ic, jc, nb, mi, ni, nq, nw, iws; + logical left; + extern logical lsame_(char *, char *); + integer nbmin, iinfo; + extern /* Subroutine */ int dorml2_(char *, char *, integer *, integer *, + integer *, doublereal *, integer *, doublereal *, doublereal *, + integer *, doublereal *, integer *), dlarfb_(char + *, char *, char *, char *, integer *, integer *, integer *, + doublereal *, integer *, doublereal *, integer *, doublereal *, + integer *, doublereal *, integer *), dlarft_(char *, char *, integer *, integer *, doublereal + *, integer *, doublereal *, doublereal *, integer *), xerbla_(char *, integer *); + extern integer ilaenv_(integer *, char *, char *, integer *, integer *, + integer *, integer *); + logical notran; + integer ldwork; + char transt[1]; + integer lwkopt; + logical lquery; + + +/* -- LAPACK routine (version 3.2) -- */ +/* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ +/* November 2006 */ + +/* .. Scalar Arguments .. */ +/* .. */ +/* .. Array Arguments .. */ +/* .. */ + +/* Purpose */ +/* ======= */ + +/* DORMLQ overwrites the general real M-by-N matrix C with */ + +/* SIDE = 'L' SIDE = 'R' */ +/* TRANS = 'N': Q * C C * Q */ +/* TRANS = 'T': Q**T * C C * Q**T */ + +/* where Q is a real orthogonal matrix defined as the product of k */ +/* elementary reflectors */ + +/* Q = H(k) . . . H(2) H(1) */ + +/* as returned by DGELQF. Q is of order M if SIDE = 'L' and of order N */ +/* if SIDE = 'R'. */ + +/* Arguments */ +/* ========= */ + +/* SIDE (input) CHARACTER*1 */ +/* = 'L': apply Q or Q**T from the Left; */ +/* = 'R': apply Q or Q**T from the Right. */ + +/* TRANS (input) CHARACTER*1 */ +/* = 'N': No transpose, apply Q; */ +/* = 'T': Transpose, apply Q**T. */ + +/* M (input) INTEGER */ +/* The number of rows of the matrix C. M >= 0. */ + +/* N (input) INTEGER */ +/* The number of columns of the matrix C. N >= 0. */ + +/* K (input) INTEGER */ +/* The number of elementary reflectors whose product defines */ +/* the matrix Q. */ +/* If SIDE = 'L', M >= K >= 0; */ +/* if SIDE = 'R', N >= K >= 0. */ + +/* A (input) DOUBLE PRECISION array, dimension */ +/* (LDA,M) if SIDE = 'L', */ +/* (LDA,N) if SIDE = 'R' */ +/* The i-th row must contain the vector which defines the */ +/* elementary reflector H(i), for i = 1,2,...,k, as returned by */ +/* DGELQF in the first k rows of its array argument A. */ +/* A is modified by the routine but restored on exit. */ + +/* LDA (input) INTEGER */ +/* The leading dimension of the array A. LDA >= max(1,K). */ + +/* TAU (input) DOUBLE PRECISION array, dimension (K) */ +/* TAU(i) must contain the scalar factor of the elementary */ +/* reflector H(i), as returned by DGELQF. */ + +/* C (input/output) DOUBLE PRECISION array, dimension (LDC,N) */ +/* On entry, the M-by-N matrix C. */ +/* On exit, C is overwritten by Q*C or Q**T*C or C*Q**T or C*Q. */ + +/* LDC (input) INTEGER */ +/* The leading dimension of the array C. LDC >= max(1,M). */ + +/* WORK (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LWORK)) */ +/* On exit, if INFO = 0, WORK(1) returns the optimal LWORK. */ + +/* LWORK (input) INTEGER */ +/* The dimension of the array WORK. */ +/* If SIDE = 'L', LWORK >= max(1,N); */ +/* if SIDE = 'R', LWORK >= max(1,M). */ +/* For optimum performance LWORK >= N*NB if SIDE = 'L', and */ +/* LWORK >= M*NB if SIDE = 'R', where NB is the optimal */ +/* blocksize. */ + +/* If LWORK = -1, then a workspace query is assumed; the routine */ +/* only calculates the optimal size of the WORK array, returns */ +/* this value as the first entry of the WORK array, and no error */ +/* message related to LWORK is issued by XERBLA. */ + +/* INFO (output) INTEGER */ +/* = 0: successful exit */ +/* < 0: if INFO = -i, the i-th argument had an illegal value */ + +/* ===================================================================== */ + +/* .. Parameters .. */ +/* .. */ +/* .. Local Scalars .. */ +/* .. */ +/* .. Local Arrays .. */ +/* .. */ +/* .. External Functions .. */ +/* .. */ +/* .. External Subroutines .. */ +/* .. */ +/* .. Intrinsic Functions .. */ +/* .. */ +/* .. Executable Statements .. */ + +/* Test the input arguments */ + + /* Parameter adjustments */ + a_dim1 = *lda; + a_offset = 1 + a_dim1; + a -= a_offset; + --tau; + c_dim1 = *ldc; + c_offset = 1 + c_dim1; + c__ -= c_offset; + --work; + + /* Function Body */ + *info = 0; + left = lsame_(side, "L"); + notran = lsame_(trans, "N"); + lquery = *lwork == -1; + +/* NQ is the order of Q and NW is the minimum dimension of WORK */ + + if (left) { + nq = *m; + nw = *n; + } else { + nq = *n; + nw = *m; + } + if (! left && ! lsame_(side, "R")) { + *info = -1; + } else if (! notran && ! lsame_(trans, "T")) { + *info = -2; + } else if (*m < 0) { + *info = -3; + } else if (*n < 0) { + *info = -4; + } else if (*k < 0 || *k > nq) { + *info = -5; + } else if (*lda < max(1,*k)) { + *info = -7; + } else if (*ldc < max(1,*m)) { + *info = -10; + } else if (*lwork < max(1,nw) && ! lquery) { + *info = -12; + } + + if (*info == 0) { + +/* Determine the block size. NB may be at most NBMAX, where NBMAX */ +/* is used to define the local array T. */ + +/* Computing MIN */ +/* Writing concatenation */ + i__3[0] = 1, a__1[0] = side; + i__3[1] = 1, a__1[1] = trans; + s_cat(ch__1, a__1, i__3, &c__2, (ftnlen)2); + i__1 = 64, i__2 = ilaenv_(&c__1, "DORMLQ", ch__1, m, n, k, &c_n1); + nb = min(i__1,i__2); + lwkopt = max(1,nw) * nb; + work[1] = (doublereal) lwkopt; + } + + if (*info != 0) { + i__1 = -(*info); + xerbla_("DORMLQ", &i__1); + return 0; + } else if (lquery) { + return 0; + } + +/* Quick return if possible */ + + if (*m == 0 || *n == 0 || *k == 0) { + work[1] = 1.; + return 0; + } + + nbmin = 2; + ldwork = nw; + if (nb > 1 && nb < *k) { + iws = nw * nb; + if (*lwork < iws) { + nb = *lwork / ldwork; +/* Computing MAX */ +/* Writing concatenation */ + i__3[0] = 1, a__1[0] = side; + i__3[1] = 1, a__1[1] = trans; + s_cat(ch__1, a__1, i__3, &c__2, (ftnlen)2); + i__1 = 2, i__2 = ilaenv_(&c__2, "DORMLQ", ch__1, m, n, k, &c_n1); + nbmin = max(i__1,i__2); + } + } else { + iws = nw; + } + + if (nb < nbmin || nb >= *k) { + +/* Use unblocked code */ + + dorml2_(side, trans, m, n, k, &a[a_offset], lda, &tau[1], &c__[ + c_offset], ldc, &work[1], &iinfo); + } else { + +/* Use blocked code */ + + if (left && notran || ! left && ! notran) { + i1 = 1; + i2 = *k; + i3 = nb; + } else { + i1 = (*k - 1) / nb * nb + 1; + i2 = 1; + i3 = -nb; + } + + if (left) { + ni = *n; + jc = 1; + } else { + mi = *m; + ic = 1; + } + + if (notran) { + *(unsigned char *)transt = 'T'; + } else { + *(unsigned char *)transt = 'N'; + } + + i__1 = i2; + i__2 = i3; + for (i__ = i1; i__2 < 0 ? i__ >= i__1 : i__ <= i__1; i__ += i__2) { +/* Computing MIN */ + i__4 = nb, i__5 = *k - i__ + 1; + ib = min(i__4,i__5); + +/* Form the triangular factor of the block reflector */ +/* H = H(i) H(i+1) . . . H(i+ib-1) */ + + i__4 = nq - i__ + 1; + dlarft_("Forward", "Rowwise", &i__4, &ib, &a[i__ + i__ * a_dim1], + lda, &tau[i__], t, &c__65); + if (left) { + +/* H or H' is applied to C(i:m,1:n) */ + + mi = *m - i__ + 1; + ic = i__; + } else { + +/* H or H' is applied to C(1:m,i:n) */ + + ni = *n - i__ + 1; + jc = i__; + } + +/* Apply H or H' */ + + dlarfb_(side, transt, "Forward", "Rowwise", &mi, &ni, &ib, &a[i__ + + i__ * a_dim1], lda, t, &c__65, &c__[ic + jc * c_dim1], + ldc, &work[1], &ldwork); +/* L10: */ + } + } + work[1] = (doublereal) lwkopt; + return 0; + +/* End of DORMLQ */ + +} /* dormlq_ */ |