LAPACK 3.3.0

cglmts.f

Go to the documentation of this file.
00001       SUBROUTINE CGLMTS( N, M, P, A, AF, LDA, B, BF, LDB, D, DF,
00002      $                   X, U, WORK, LWORK, RWORK, RESULT )
00003 *
00004 *  -- LAPACK test routine (version 3.1) --
00005 *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
00006 *     November 2006
00007 *
00008 *     .. Scalar Arguments ..
00009       INTEGER            LDA, LDB, LWORK, M, P, N
00010       REAL               RESULT
00011 *     ..
00012 *     .. Array Arguments ..
00013       REAL               RWORK( * )
00014       COMPLEX            A( LDA, * ), AF( LDA, * ), B( LDB, * ),
00015      $                   BF( LDB, * ), D( * ), DF( * ), U( * ),
00016      $                   WORK( LWORK ), X( * )
00017 *
00018 *  Purpose
00019 *  =======
00020 *
00021 *  CGLMTS tests CGGGLM - a subroutine for solving the generalized
00022 *  linear model problem.
00023 *
00024 *  Arguments
00025 *  =========
00026 *
00027 *  N       (input) INTEGER
00028 *          The number of rows of the matrices A and B.  N >= 0.
00029 *
00030 *  M       (input) INTEGER
00031 *          The number of columns of the matrix A.  M >= 0.
00032 *
00033 *  P       (input) INTEGER
00034 *          The number of columns of the matrix B.  P >= 0.
00035 *
00036 *  A       (input) COMPLEX array, dimension (LDA,M)
00037 *          The N-by-M matrix A.
00038 *
00039 *  AF      (workspace) COMPLEX array, dimension (LDA,M)
00040 *
00041 *  LDA     (input) INTEGER
00042 *          The leading dimension of the arrays A, AF. LDA >= max(M,N).
00043 *
00044 *  B       (input) COMPLEX array, dimension (LDB,P)
00045 *          The N-by-P matrix A.
00046 *
00047 *  BF      (workspace) COMPLEX array, dimension (LDB,P)
00048 *
00049 *  LDB     (input) INTEGER
00050 *          The leading dimension of the arrays B, BF. LDB >= max(P,N).
00051 *
00052 *  D       (input) COMPLEX array, dimension( N )
00053 *          On input, the left hand side of the GLM.
00054 *
00055 *  DF      (workspace) COMPLEX array, dimension( N )
00056 *
00057 *  X       (output) COMPLEX array, dimension( M )
00058 *          solution vector X in the GLM problem.
00059 *
00060 *  U       (output) COMPLEX array, dimension( P )
00061 *          solution vector U in the GLM problem.
00062 *
00063 *  WORK    (workspace) COMPLEX array, dimension (LWORK)
00064 *
00065 *  LWORK   (input) INTEGER
00066 *          The dimension of the array WORK.
00067 *
00068 *  RWORK   (workspace) REAL array, dimension (M)
00069 *
00070 *  RESULT   (output) REAL
00071 *          The test ratio:
00072 *                           norm( d - A*x - B*u )
00073 *            RESULT = -----------------------------------------
00074 *                     (norm(A)+norm(B))*(norm(x)+norm(u))*EPS
00075 *
00076 *  ====================================================================
00077 *
00078 *     .. Parameters ..
00079       REAL               ZERO
00080       PARAMETER          ( ZERO = 0.0E+0 )
00081       COMPLEX            CONE
00082       PARAMETER          ( CONE = 1.0E+0 )
00083 *     ..
00084 *     .. Local Scalars ..
00085       INTEGER            INFO
00086       REAL               ANORM, BNORM, EPS, XNORM, YNORM, DNORM, UNFL
00087 *     ..
00088 *     .. External Functions ..
00089       REAL               SCASUM, SLAMCH, CLANGE
00090       EXTERNAL           SCASUM, SLAMCH, CLANGE
00091 *     ..
00092 *     .. External Subroutines ..
00093       EXTERNAL           CLACPY
00094 *
00095 *     .. Intrinsic Functions ..
00096       INTRINSIC          MAX
00097 *     ..
00098 *     .. Executable Statements ..
00099 *
00100       EPS = SLAMCH( 'Epsilon' )
00101       UNFL = SLAMCH( 'Safe minimum' )
00102       ANORM = MAX( CLANGE( '1', N, M, A, LDA, RWORK ), UNFL )
00103       BNORM = MAX( CLANGE( '1', N, P, B, LDB, RWORK ), UNFL )
00104 *
00105 *     Copy the matrices A and B to the arrays AF and BF,
00106 *     and the vector D the array DF.
00107 *
00108       CALL CLACPY( 'Full', N, M, A, LDA, AF, LDA )
00109       CALL CLACPY( 'Full', N, P, B, LDB, BF, LDB )
00110       CALL CCOPY( N, D, 1, DF, 1 )
00111 *
00112 *     Solve GLM problem
00113 *
00114       CALL CGGGLM( N, M, P, AF, LDA, BF, LDB, DF, X, U, WORK, LWORK,
00115      $             INFO )
00116 *
00117 *     Test the residual for the solution of LSE
00118 *
00119 *                       norm( d - A*x - B*u )
00120 *       RESULT = -----------------------------------------
00121 *                (norm(A)+norm(B))*(norm(x)+norm(u))*EPS
00122 *
00123       CALL CCOPY( N, D, 1, DF, 1 )
00124       CALL CGEMV( 'No transpose', N, M, -CONE, A, LDA, X, 1, CONE,
00125      $            DF, 1 )
00126 *
00127       CALL CGEMV( 'No transpose', N, P, -CONE, B, LDB, U, 1, CONE,
00128      $            DF, 1 )
00129 *
00130       DNORM = SCASUM( N, DF, 1 )
00131       XNORM = SCASUM( M, X, 1 ) + SCASUM( P, U, 1 )
00132       YNORM = ANORM + BNORM
00133 *
00134       IF( XNORM.LE.ZERO ) THEN
00135          RESULT = ZERO
00136       ELSE
00137          RESULT =  ( ( DNORM / YNORM ) / XNORM ) /EPS
00138       END IF
00139 *
00140       RETURN
00141 *
00142 *     End of CGLMTS
00143 *
00144       END
 All Files Functions