Files
lapack/TESTING/LIN/dlqt01.f
Julie e18d437924 Making LAPACK's code eternal... no more version and date in source files.
GitHub is now enabling us to track accurately version and date.
No need for this anymore.
2021-03-25 10:16:58 -07:00

230 lines
6.0 KiB
FortranFixed

*> \brief \b DLQT01
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
* Definition:
* ===========
*
* SUBROUTINE DLQT01( M, N, A, AF, Q, L, LDA, TAU, WORK, LWORK,
* RWORK, RESULT )
*
* .. Scalar Arguments ..
* INTEGER LDA, LWORK, M, N
* ..
* .. Array Arguments ..
* DOUBLE PRECISION A( LDA, * ), AF( LDA, * ), L( LDA, * ),
* $ Q( LDA, * ), RESULT( * ), RWORK( * ), TAU( * ),
* $ WORK( LWORK )
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DLQT01 tests DGELQF, which computes the LQ factorization of an m-by-n
*> matrix A, and partially tests DORGLQ which forms the n-by-n
*> orthogonal matrix Q.
*>
*> DLQT01 compares L with A*Q', and checks that Q is orthogonal.
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] M
*> \verbatim
*> M is INTEGER
*> The number of rows of the matrix A. M >= 0.
*> \endverbatim
*>
*> \param[in] N
*> \verbatim
*> N is INTEGER
*> The number of columns of the matrix A. N >= 0.
*> \endverbatim
*>
*> \param[in] A
*> \verbatim
*> A is DOUBLE PRECISION array, dimension (LDA,N)
*> The m-by-n matrix A.
*> \endverbatim
*>
*> \param[out] AF
*> \verbatim
*> AF is DOUBLE PRECISION array, dimension (LDA,N)
*> Details of the LQ factorization of A, as returned by DGELQF.
*> See DGELQF for further details.
*> \endverbatim
*>
*> \param[out] Q
*> \verbatim
*> Q is DOUBLE PRECISION array, dimension (LDA,N)
*> The n-by-n orthogonal matrix Q.
*> \endverbatim
*>
*> \param[out] L
*> \verbatim
*> L is DOUBLE PRECISION array, dimension (LDA,max(M,N))
*> \endverbatim
*>
*> \param[in] LDA
*> \verbatim
*> LDA is INTEGER
*> The leading dimension of the arrays A, AF, Q and L.
*> LDA >= max(M,N).
*> \endverbatim
*>
*> \param[out] TAU
*> \verbatim
*> TAU is DOUBLE PRECISION array, dimension (min(M,N))
*> The scalar factors of the elementary reflectors, as returned
*> by DGELQF.
*> \endverbatim
*>
*> \param[out] WORK
*> \verbatim
*> WORK is DOUBLE PRECISION array, dimension (LWORK)
*> \endverbatim
*>
*> \param[in] LWORK
*> \verbatim
*> LWORK is INTEGER
*> The dimension of the array WORK.
*> \endverbatim
*>
*> \param[out] RWORK
*> \verbatim
*> RWORK is DOUBLE PRECISION array, dimension (max(M,N))
*> \endverbatim
*>
*> \param[out] RESULT
*> \verbatim
*> RESULT is DOUBLE PRECISION array, dimension (2)
*> The test ratios:
*> RESULT(1) = norm( L - A*Q' ) / ( N * norm(A) * EPS )
*> RESULT(2) = norm( I - Q*Q' ) / ( N * EPS )
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \ingroup double_lin
*
* =====================================================================
SUBROUTINE DLQT01( M, N, A, AF, Q, L, LDA, TAU, WORK, LWORK,
$ RWORK, RESULT )
*
* -- LAPACK test routine --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
* .. Scalar Arguments ..
INTEGER LDA, LWORK, M, N
* ..
* .. Array Arguments ..
DOUBLE PRECISION A( LDA, * ), AF( LDA, * ), L( LDA, * ),
$ Q( LDA, * ), RESULT( * ), RWORK( * ), TAU( * ),
$ WORK( LWORK )
* ..
*
* =====================================================================
*
* .. Parameters ..
DOUBLE PRECISION ZERO, ONE
PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
DOUBLE PRECISION ROGUE
PARAMETER ( ROGUE = -1.0D+10 )
* ..
* .. Local Scalars ..
INTEGER INFO, MINMN
DOUBLE PRECISION ANORM, EPS, RESID
* ..
* .. External Functions ..
DOUBLE PRECISION DLAMCH, DLANGE, DLANSY
EXTERNAL DLAMCH, DLANGE, DLANSY
* ..
* .. External Subroutines ..
EXTERNAL DGELQF, DGEMM, DLACPY, DLASET, DORGLQ, DSYRK
* ..
* .. Intrinsic Functions ..
INTRINSIC DBLE, MAX, MIN
* ..
* .. Scalars in Common ..
CHARACTER*32 SRNAMT
* ..
* .. Common blocks ..
COMMON / SRNAMC / SRNAMT
* ..
* .. Executable Statements ..
*
MINMN = MIN( M, N )
EPS = DLAMCH( 'Epsilon' )
*
* Copy the matrix A to the array AF.
*
CALL DLACPY( 'Full', M, N, A, LDA, AF, LDA )
*
* Factorize the matrix A in the array AF.
*
SRNAMT = 'DGELQF'
CALL DGELQF( M, N, AF, LDA, TAU, WORK, LWORK, INFO )
*
* Copy details of Q
*
CALL DLASET( 'Full', N, N, ROGUE, ROGUE, Q, LDA )
IF( N.GT.1 )
$ CALL DLACPY( 'Upper', M, N-1, AF( 1, 2 ), LDA, Q( 1, 2 ), LDA )
*
* Generate the n-by-n matrix Q
*
SRNAMT = 'DORGLQ'
CALL DORGLQ( N, N, MINMN, Q, LDA, TAU, WORK, LWORK, INFO )
*
* Copy L
*
CALL DLASET( 'Full', M, N, ZERO, ZERO, L, LDA )
CALL DLACPY( 'Lower', M, N, AF, LDA, L, LDA )
*
* Compute L - A*Q'
*
CALL DGEMM( 'No transpose', 'Transpose', M, N, N, -ONE, A, LDA, Q,
$ LDA, ONE, L, LDA )
*
* Compute norm( L - Q'*A ) / ( N * norm(A) * EPS ) .
*
ANORM = DLANGE( '1', M, N, A, LDA, RWORK )
RESID = DLANGE( '1', M, N, L, LDA, RWORK )
IF( ANORM.GT.ZERO ) THEN
RESULT( 1 ) = ( ( RESID / DBLE( MAX( 1, N ) ) ) / ANORM ) / EPS
ELSE
RESULT( 1 ) = ZERO
END IF
*
* Compute I - Q*Q'
*
CALL DLASET( 'Full', N, N, ZERO, ONE, L, LDA )
CALL DSYRK( 'Upper', 'No transpose', N, N, -ONE, Q, LDA, ONE, L,
$ LDA )
*
* Compute norm( I - Q*Q' ) / ( N * EPS ) .
*
RESID = DLANSY( '1', 'Upper', N, L, LDA, RWORK )
*
RESULT( 2 ) = ( RESID / DBLE( MAX( 1, N ) ) ) / EPS
*
RETURN
*
* End of DLQT01
*
END