TESTING/EIG/dlsets.f

   1 *> \brief \b DLSETS
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *  Definition:
   9 *  ===========
  10 *
  11 *       SUBROUTINE DLSETS( M, P, N, A, AF, LDA, B, BF, LDB, C, CF, D, DF,
  12 *                          X, WORK, LWORK, RWORK, RESULT )
  13 *
  14 *       .. Scalar Arguments ..
  15 *       INTEGER            LDA, LDB, LWORK, M, N, P
  16 *       ..
  17 *       .. Array Arguments ..
  18 *
  19 *
  20 *> \par Purpose:
  21 *  =============
  22 *>
  23 *> \verbatim
  24 *>
  25 *> DLSETS tests DGGLSE - a subroutine for solving linear equality
  26 *> constrained least square problem (LSE).
  27 *> \endverbatim
  28 *
  29 *  Arguments:
  30 *  ==========
  31 *
  32 *> \param[in] M
  33 *> \verbatim
  34 *>          M is INTEGER
  35 *>          The number of rows of the matrix A.  M >= 0.
  36 *> \endverbatim
  37 *>
  38 *> \param[in] P
  39 *> \verbatim
  40 *>          P is INTEGER
  41 *>          The number of rows of the matrix B.  P >= 0.
  42 *> \endverbatim
  43 *>
  44 *> \param[in] N
  45 *> \verbatim
  46 *>          N is INTEGER
  47 *>          The number of columns of the matrices A and B.  N >= 0.
  48 *> \endverbatim
  49 *>
  50 *> \param[in] A
  51 *> \verbatim
  52 *>          A is DOUBLE PRECISION array, dimension (LDA,N)
  53 *>          The M-by-N matrix A.
  54 *> \endverbatim
  55 *>
  56 *> \param[out] AF
  57 *> \verbatim
  58 *>          AF is DOUBLE PRECISION array, dimension (LDA,N)
  59 *> \endverbatim
  60 *>
  61 *> \param[in] LDA
  62 *> \verbatim
  63 *>          LDA is INTEGER
  64 *>          The leading dimension of the arrays A, AF, Q and R.
  65 *>          LDA >= max(M,N).
  66 *> \endverbatim
  67 *>
  68 *> \param[in] B
  69 *> \verbatim
  70 *>          B is DOUBLE PRECISION array, dimension (LDB,N)
  71 *>          The P-by-N matrix A.
  72 *> \endverbatim
  73 *>
  74 *> \param[out] BF
  75 *> \verbatim
  76 *>          BF is DOUBLE PRECISION array, dimension (LDB,N)
  77 *> \endverbatim
  78 *>
  79 *> \param[in] LDB
  80 *> \verbatim
  81 *>          LDB is INTEGER
  82 *>          The leading dimension of the arrays B, BF, V and S.
  83 *>          LDB >= max(P,N).
  84 *> \endverbatim
  85 *>
  86 *> \param[in] C
  87 *> \verbatim
  88 *>          C is DOUBLE PRECISION array, dimension( M )
  89 *>          the vector C in the LSE problem.
  90 *> \endverbatim
  91 *>
  92 *> \param[out] CF
  93 *> \verbatim
  94 *>          CF is DOUBLE PRECISION array, dimension( M )
  95 *> \endverbatim
  96 *>
  97 *> \param[in] D
  98 *> \verbatim
  99 *>          D is DOUBLE PRECISION array, dimension( P )
 100 *>          the vector D in the LSE problem.
 101 *> \endverbatim
 102 *>
 103 *> \param[out] DF
 104 *> \verbatim
 105 *>          DF is DOUBLE PRECISION array, dimension( P )
 106 *> \endverbatim
 107 *>
 108 *> \param[out] X
 109 *> \verbatim
 110 *>          X is DOUBLE PRECISION array, dimension( N )
 111 *>          solution vector X in the LSE problem.
 112 *> \endverbatim
 113 *>
 114 *> \param[out] WORK
 115 *> \verbatim
 116 *>          WORK is DOUBLE PRECISION array, dimension (LWORK)
 117 *> \endverbatim
 118 *>
 119 *> \param[in] LWORK
 120 *> \verbatim
 121 *>          LWORK is INTEGER
 122 *>          The dimension of the array WORK.
 123 *> \endverbatim
 124 *>
 125 *> \param[out] RWORK
 126 *> \verbatim
 127 *>          RWORK is DOUBLE PRECISION array, dimension (M)
 128 *> \endverbatim
 129 *>
 130 *> \param[out] RESULT
 131 *> \verbatim
 132 *>          RESULT is DOUBLE PRECISION array, dimension (2)
 133 *>          The test ratios:
 134 *>            RESULT(1) = norm( A*x - c )/ norm(A)*norm(X)*EPS
 135 *>            RESULT(2) = norm( B*x - d )/ norm(B)*norm(X)*EPS
 136 *> \endverbatim
 137 *
 138 *  Authors:
 139 *  ========
 140 *
 141 *> \author Univ. of Tennessee
 142 *> \author Univ. of California Berkeley
 143 *> \author Univ. of Colorado Denver
 144 *> \author NAG Ltd.
 145 *
 146 *> \date November 2011
 147 *
 148 *> \ingroup double_eig
 149 *
 150 *  =====================================================================
 151       SUBROUTINE DLSETS( M, P, N, A, AF, LDA, B, BF, LDB, C, CF, D, DF,
 152      $                   X, WORK, LWORK, RWORK, RESULT )
 153 *
 154 *  -- LAPACK test routine (version 3.4.0) --
 155 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 156 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 157 *     November 2011
 158 *
 159 *     .. Scalar Arguments ..
 160       INTEGER            LDA, LDB, LWORK, M, N, P
 161 *     ..
 162 *     .. Array Arguments ..
 163 *
 164 *  ====================================================================
 165 *
 166       DOUBLE PRECISION   A( LDA, * ), AF( LDA, * ), B( LDB, * ),
 167      $                   BF( LDB, * ), C( * ), CF( * ), D( * ), DF( * ),
 168      $                   RESULT( 2 ), RWORK( * ), WORK( LWORK ), X( * )
 169 *     ..
 170 *     .. Local Scalars ..
 171       INTEGER            INFO
 172 *     ..
 173 *     .. External Subroutines ..
 174       EXTERNAL           DCOPY, DGET02, DGGLSE, DLACPY
 175 *     ..
 176 *     .. Executable Statements ..
 177 *
 178 *     Copy the matrices A and B to the arrays AF and BF,
 179 *     and the vectors C and D to the arrays CF and DF,
 180 *
 181       CALL DLACPY( 'Full', M, N, A, LDA, AF, LDA )
 182       CALL DLACPY( 'Full', P, N, B, LDB, BF, LDB )
 183       CALL DCOPY( M, C, 1, CF, 1 )
 184       CALL DCOPY( P, D, 1, DF, 1 )
 185 *
 186 *     Solve LSE problem
 187 *
 188       CALL DGGLSE( M, N, P, AF, LDA, BF, LDB, CF, DF, X, WORK, LWORK,
 189      $             INFO )
 190 *
 191 *     Test the residual for the solution of LSE
 192 *
 193 *     Compute RESULT(1) = norm( A*x - c ) / norm(A)*norm(X)*EPS
 194 *
 195       CALL DCOPY( M, C, 1, CF, 1 )
 196       CALL DCOPY( P, D, 1, DF, 1 )
 197       CALL DGET02( 'No transpose', M, N, 1, A, LDA, X, N, CF, M, RWORK,
 198      $             RESULT( 1 ) )
 199 *
 200 *     Compute result(2) = norm( B*x - d ) / norm(B)*norm(X)*EPS
 201 *
 202       CALL DGET02( 'No transpose', P, N, 1, B, LDB, X, N, DF, P, RWORK,
 203      $             RESULT( 2 ) )
 204 *
 205       RETURN
 206 *
 207 *     End of DLSETS
 208 *
 209       END