LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dglmts.f
Go to the documentation of this file.
1*> \brief \b DGLMTS
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 DGLMTS( N, M, P, A, AF, LDA, B, BF, LDB, D, DF, X, U,
12* WORK, LWORK, RWORK, RESULT )
13*
14* .. Scalar Arguments ..
15* INTEGER LDA, LDB, LWORK, M, N, P
16* DOUBLE PRECISION RESULT
17* ..
18* .. Array Arguments ..
19*
20*
21*> \par Purpose:
22* =============
23*>
24*> \verbatim
25*>
26*> DGLMTS tests DGGGLM - a subroutine for solving the generalized
27*> linear model problem.
28*> \endverbatim
29*
30* Arguments:
31* ==========
32*
33*> \param[in] N
34*> \verbatim
35*> N is INTEGER
36*> The number of rows of the matrices A and B. N >= 0.
37*> \endverbatim
38*>
39*> \param[in] M
40*> \verbatim
41*> M is INTEGER
42*> The number of columns of the matrix A. M >= 0.
43*> \endverbatim
44*>
45*> \param[in] P
46*> \verbatim
47*> P is INTEGER
48*> The number of columns of the matrix B. P >= 0.
49*> \endverbatim
50*>
51*> \param[in] A
52*> \verbatim
53*> A is DOUBLE PRECISION array, dimension (LDA,M)
54*> The N-by-M matrix A.
55*> \endverbatim
56*>
57*> \param[out] AF
58*> \verbatim
59*> AF is DOUBLE PRECISION array, dimension (LDA,M)
60*> \endverbatim
61*>
62*> \param[in] LDA
63*> \verbatim
64*> LDA is INTEGER
65*> The leading dimension of the arrays A, AF. LDA >= max(M,N).
66*> \endverbatim
67*>
68*> \param[in] B
69*> \verbatim
70*> B is DOUBLE PRECISION array, dimension (LDB,P)
71*> The N-by-P matrix A.
72*> \endverbatim
73*>
74*> \param[out] BF
75*> \verbatim
76*> BF is DOUBLE PRECISION array, dimension (LDB,P)
77*> \endverbatim
78*>
79*> \param[in] LDB
80*> \verbatim
81*> LDB is INTEGER
82*> The leading dimension of the arrays B, BF. LDB >= max(P,N).
83*> \endverbatim
84*>
85*> \param[in] D
86*> \verbatim
87*> D is DOUBLE PRECISION array, dimension( N )
88*> On input, the left hand side of the GLM.
89*> \endverbatim
90*>
91*> \param[out] DF
92*> \verbatim
93*> DF is DOUBLE PRECISION array, dimension( N )
94*> \endverbatim
95*>
96*> \param[out] X
97*> \verbatim
98*> X is DOUBLE PRECISION array, dimension( M )
99*> solution vector X in the GLM problem.
100*> \endverbatim
101*>
102*> \param[out] U
103*> \verbatim
104*> U is DOUBLE PRECISION array, dimension( P )
105*> solution vector U in the GLM problem.
106*> \endverbatim
107*>
108*> \param[out] WORK
109*> \verbatim
110*> WORK is DOUBLE PRECISION array, dimension (LWORK)
111*> \endverbatim
112*>
113*> \param[in] LWORK
114*> \verbatim
115*> LWORK is INTEGER
116*> The dimension of the array WORK.
117*> \endverbatim
118*>
119*> \param[out] RWORK
120*> \verbatim
121*> RWORK is DOUBLE PRECISION array, dimension (M)
122*> \endverbatim
123*>
124*> \param[out] RESULT
125*> \verbatim
126*> RESULT is DOUBLE PRECISION
127*> The test ratio:
128*> norm( d - A*x - B*u )
129*> RESULT = -----------------------------------------
130*> (norm(A)+norm(B))*(norm(x)+norm(u))*EPS
131*> \endverbatim
132*
133* Authors:
134* ========
135*
136*> \author Univ. of Tennessee
137*> \author Univ. of California Berkeley
138*> \author Univ. of Colorado Denver
139*> \author NAG Ltd.
140*
141*> \ingroup double_eig
142*
143* =====================================================================
144 SUBROUTINE dglmts( N, M, P, A, AF, LDA, B, BF, LDB, D, DF, X, U,
145 $ WORK, LWORK, RWORK, RESULT )
146*
147* -- LAPACK test routine --
148* -- LAPACK is a software package provided by Univ. of Tennessee, --
149* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
150*
151* .. Scalar Arguments ..
152 INTEGER LDA, LDB, LWORK, M, N, P
153 DOUBLE PRECISION RESULT
154* ..
155* .. Array Arguments ..
156*
157* ====================================================================
158*
159 DOUBLE PRECISION A( LDA, * ), AF( LDA, * ), B( LDB, * ),
160 $ bf( ldb, * ), d( * ), df( * ), rwork( * ),
161 $ u( * ), work( lwork ), x( * )
162* ..
163* .. Parameters ..
164 DOUBLE PRECISION ZERO, ONE
165 parameter( zero = 0.0d+0, one = 1.0d+0 )
166* ..
167* .. Local Scalars ..
168 INTEGER INFO
169 DOUBLE PRECISION ANORM, BNORM, DNORM, EPS, UNFL, XNORM, YNORM
170* ..
171* .. External Functions ..
172 DOUBLE PRECISION DASUM, DLAMCH, DLANGE
173 EXTERNAL dasum, dlamch, dlange
174* ..
175* .. External Subroutines ..
176*
177 EXTERNAL dcopy, dgemv, dggglm, dlacpy
178* ..
179* .. Intrinsic Functions ..
180 INTRINSIC max
181* ..
182* .. Executable Statements ..
183*
184 eps = dlamch( 'Epsilon' )
185 unfl = dlamch( 'Safe minimum' )
186 anorm = max( dlange( '1', n, m, a, lda, rwork ), unfl )
187 bnorm = max( dlange( '1', n, p, b, ldb, rwork ), unfl )
188*
189* Copy the matrices A and B to the arrays AF and BF,
190* and the vector D the array DF.
191*
192 CALL dlacpy( 'Full', n, m, a, lda, af, lda )
193 CALL dlacpy( 'Full', n, p, b, ldb, bf, ldb )
194 CALL dcopy( n, d, 1, df, 1 )
195*
196* Solve GLM problem
197*
198 CALL dggglm( n, m, p, af, lda, bf, ldb, df, x, u, work, lwork,
199 $ info )
200*
201* Test the residual for the solution of LSE
202*
203* norm( d - A*x - B*u )
204* RESULT = -----------------------------------------
205* (norm(A)+norm(B))*(norm(x)+norm(u))*EPS
206*
207 CALL dcopy( n, d, 1, df, 1 )
208 CALL dgemv( 'No transpose', n, m, -one, a, lda, x, 1, one, df, 1 )
209*
210 CALL dgemv( 'No transpose', n, p, -one, b, ldb, u, 1, one, df, 1 )
211*
212 dnorm = dasum( n, df, 1 )
213 xnorm = dasum( m, x, 1 ) + dasum( p, u, 1 )
214 ynorm = anorm + bnorm
215*
216 IF( xnorm.LE.zero ) THEN
217 result = zero
218 ELSE
219 result = ( ( dnorm / ynorm ) / xnorm ) / eps
220 END IF
221*
222 RETURN
223*
224* End of DGLMTS
225*
226 END
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:82
subroutine dgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DGEMV
Definition: dgemv.f:156
subroutine dglmts(N, M, P, A, AF, LDA, B, BF, LDB, D, DF, X, U, WORK, LWORK, RWORK, RESULT)
DGLMTS
Definition: dglmts.f:146
subroutine dggglm(N, M, P, A, LDA, B, LDB, D, X, Y, WORK, LWORK, INFO)
DGGGLM
Definition: dggglm.f:185