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

   1 #include "clapack.h"
   2
   3 /* Table of constant values */
   4
   5 static doublereal c_b9 = 1.;
   6
   7 /* Subroutine */ int dpotrs_(char *uplo, integer *n, integer *nrhs,
   8         doublereal *a, integer *lda, doublereal *b, integer *ldb, integer *
   9         info)
  10 {
  11     /* System generated locals */
  12     integer a_dim1, a_offset, b_dim1, b_offset, i__1;
  13
  14     /* Local variables */
  15     extern logical lsame_(char *, char *);
  16     extern /* Subroutine */ int dtrsm_(char *, char *, char *, char *,
  17             integer *, integer *, doublereal *, doublereal *, integer *,
  18             doublereal *, integer *);
  19     logical upper;
  20     extern /* Subroutine */ int xerbla_(char *, integer *);
  21
  22
  23 /*  -- LAPACK routine (version 3.1) -- */
  24 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
  25 /*     November 2006 */
  26
  27 /*     .. Scalar Arguments .. */
  28 /*     .. */
  29 /*     .. Array Arguments .. */
  30 /*     .. */
  31
  32 /*  Purpose */
  33 /*  ======= */
  34
  35 /*  DPOTRS solves a system of linear equations A*X = B with a symmetric */
  36 /*  positive definite matrix A using the Cholesky factorization */
  37 /*  A = U**T*U or A = L*L**T computed by DPOTRF. */
  38
  39 /*  Arguments */
  40 /*  ========= */
  41
  42 /*  UPLO    (input) CHARACTER*1 */
  43 /*          = 'U':  Upper triangle of A is stored; */
  44 /*          = 'L':  Lower triangle of A is stored. */
  45
  46 /*  N       (input) INTEGER */
  47 /*          The order of the matrix A.  N >= 0. */
  48
  49 /*  NRHS    (input) INTEGER */
  50 /*          The number of right hand sides, i.e., the number of columns */
  51 /*          of the matrix B.  NRHS >= 0. */
  52
  53 /*  A       (input) DOUBLE PRECISION array, dimension (LDA,N) */
  54 /*          The triangular factor U or L from the Cholesky factorization */
  55 /*          A = U**T*U or A = L*L**T, as computed by DPOTRF. */
  56
  57 /*  LDA     (input) INTEGER */
  58 /*          The leading dimension of the array A.  LDA >= max(1,N). */
  59
  60 /*  B       (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS) */
  61 /*          On entry, the right hand side matrix B. */
  62 /*          On exit, the solution matrix X. */
  63
  64 /*  LDB     (input) INTEGER */
  65 /*          The leading dimension of the array B.  LDB >= max(1,N). */
  66
  67 /*  INFO    (output) INTEGER */
  68 /*          = 0:  successful exit */
  69 /*          < 0:  if INFO = -i, the i-th argument had an illegal value */
  70
  71 /*  ===================================================================== */
  72
  73 /*     .. Parameters .. */
  74 /*     .. */
  75 /*     .. Local Scalars .. */
  76 /*     .. */
  77 /*     .. External Functions .. */
  78 /*     .. */
  79 /*     .. External Subroutines .. */
  80 /*     .. */
  81 /*     .. Intrinsic Functions .. */
  82 /*     .. */
  83 /*     .. Executable Statements .. */
  84
  85 /*     Test the input parameters. */
  86
  87     /* Parameter adjustments */
  88     a_dim1 = *lda;
  89     a_offset = 1 + a_dim1;
  90     a -= a_offset;
  91     b_dim1 = *ldb;
  92     b_offset = 1 + b_dim1;
  93     b -= b_offset;
  94
  95     /* Function Body */
  96     *info = 0;
  97     upper = lsame_(uplo, "U");
  98     if (! upper && ! lsame_(uplo, "L")) {
  99         *info = -1;
 100     } else if (*n < 0) {
 101         *info = -2;
 102     } else if (*nrhs < 0) {
 103         *info = -3;
 104     } else if (*lda < max(1,*n)) {
 105         *info = -5;
 106     } else if (*ldb < max(1,*n)) {
 107         *info = -7;
 108     }
 109     if (*info != 0) {
 110         i__1 = -(*info);
 111         xerbla_("DPOTRS", &i__1);
 112         return 0;
 113     }
 114
 115 /*     Quick return if possible */
 116
 117     if (*n == 0 || *nrhs == 0) {
 118         return 0;
 119     }
 120
 121     if (upper) {
 122
 123 /*        Solve A*X = B where A = U'*U. */
 124
 125 /*        Solve U'*X = B, overwriting B with X. */
 126
 127         dtrsm_("Left", "Upper", "Transpose", "Non-unit", n, nrhs, &c_b9, &a[
 128                 a_offset], lda, &b[b_offset], ldb);
 129
 130 /*        Solve U*X = B, overwriting B with X. */
 131
 132         dtrsm_("Left", "Upper", "No transpose", "Non-unit", n, nrhs, &c_b9, &
 133                 a[a_offset], lda, &b[b_offset], ldb);
 134     } else {
 135
 136 /*        Solve A*X = B where A = L*L'. */
 137
 138 /*        Solve L*X = B, overwriting B with X. */
 139
 140         dtrsm_("Left", "Lower", "No transpose", "Non-unit", n, nrhs, &c_b9, &
 141                 a[a_offset], lda, &b[b_offset], ldb);
 142
 143 /*        Solve L'*X = B, overwriting B with X. */
 144
 145         dtrsm_("Left", "Lower", "Transpose", "Non-unit", n, nrhs, &c_b9, &a[
 146                 a_offset], lda, &b[b_offset], ldb);
 147     }
 148
 149     return 0;
 150
 151 /*     End of DPOTRS */
 152
 153 } /* dpotrs_ */