LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ zsyl01()

subroutine zsyl01 ( double precision  thresh,
integer, dimension( 3 )  nfail,
double precision, dimension( 2 )  rmax,
integer, dimension( 2 )  ninfo,
integer  knt 
)

ZSYL01

Purpose:
 ZSYL01 tests ZTRSYL and ZTRSYL3, routines for solving the Sylvester matrix
 equation

    op(A)*X + ISGN*X*op(B) = scale*C,

 where op(A) and op(B) are both upper triangular form, op() represents an
 optional conjugate transpose, and ISGN can be -1 or +1. Scale is an output
 less than or equal to 1, chosen to avoid overflow in X.

 The test code verifies that the following residual does not exceed
 the provided threshold:

    norm(op(A)*X + ISGN*X*op(B) - scale*C) /
        (EPS*max(norm(A),norm(B))*norm(X))

 This routine complements ZGET35 by testing with larger,
 random matrices, of which some require rescaling of X to avoid overflow.
Parameters
[in]THRESH
          THRESH is DOUBLE PRECISION
          A test will count as "failed" if the residual, computed as
          described above, exceeds THRESH.
[out]NFAIL
          NFAIL is INTEGER array, dimension (3)
          NFAIL(1) = No. of times residual ZTRSYL exceeds threshold THRESH
          NFAIL(2) = No. of times residual ZTRSYL3 exceeds threshold THRESH
          NFAIL(3) = No. of times ZTRSYL3 and ZTRSYL deviate
[out]RMAX
          RMAX is DOUBLE PRECISION array, dimension (2)
          RMAX(1) = Value of the largest test ratio of ZTRSYL
          RMAX(2) = Value of the largest test ratio of ZTRSYL3
[out]NINFO
          NINFO is INTEGER array, dimension (2)
          NINFO(1) = No. of times ZTRSYL returns an expected INFO
          NINFO(2) = No. of times ZTRSYL3 returns an expected INFO
[out]KNT
          KNT is INTEGER
          Total number of examples tested.

Definition at line 88 of file zsyl01.f.

89 IMPLICIT NONE
90*
91* -- LAPACK test routine --
92* -- LAPACK is a software package provided by Univ. of Tennessee, --
93* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
94*
95* .. Scalar Arguments ..
96 INTEGER KNT
97 DOUBLE PRECISION THRESH
98* ..
99* .. Array Arguments ..
100 INTEGER NFAIL( 3 ), NINFO( 2 )
101 DOUBLE PRECISION RMAX( 2 )
102* ..
103*
104* =====================================================================
105* ..
106* .. Parameters ..
107 COMPLEX*16 CONE
108 parameter( cone = ( 1.0d0, 0.0d+0 ) )
109 DOUBLE PRECISION ONE, ZERO
110 parameter( zero = 0.0d+0, one = 1.0d+0 )
111 INTEGER MAXM, MAXN, LDSWORK
112 parameter( maxm = 185, maxn = 192, ldswork = 36 )
113* ..
114* .. Local Scalars ..
115 CHARACTER TRANA, TRANB
116 INTEGER I, INFO, IINFO, ISGN, ITRANA, ITRANB, J, KLA,
117 $ KUA, KLB, KUB, M, N
118 DOUBLE PRECISION ANRM, BNRM, BIGNUM, EPS, RES, RES1,
119 $ SCALE, SCALE3, SMLNUM, TNRM, XNRM
120 COMPLEX*16 RMUL
121* ..
122* .. Local Arrays ..
123 COMPLEX*16 DUML( MAXM ), DUMR( MAXN ),
124 $ D( MAX( MAXM, MAXN ) )
125 DOUBLE PRECISION DUM( MAXN ), VM( 2 )
126 INTEGER ISEED( 4 ), IWORK( MAXM + MAXN + 2 )
127* ..
128* .. Allocatable Arrays ..
129 INTEGER AllocateStatus
130 COMPLEX*16, DIMENSION(:,:), ALLOCATABLE :: A, B, C, CC, X
131 DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: SWORK
132* ..
133* .. External Functions ..
134 LOGICAL DISNAN
135 DOUBLE PRECISION DLAMCH, ZLANGE
136 EXTERNAL disnan, dlamch, zlange
137* ..
138* .. External Subroutines ..
139 EXTERNAL zlatmr, zlacpy, zgemm, ztrsyl, ztrsyl3
140* ..
141* .. Intrinsic Functions ..
142 INTRINSIC abs, dble, max, sqrt
143* ..
144* .. Allocate memory dynamically ..
145 ALLOCATE ( a( maxm, maxm ), stat = allocatestatus )
146 IF( allocatestatus /= 0 ) stop "*** Not enough memory ***"
147 ALLOCATE ( b( maxn, maxn ), stat = allocatestatus )
148 IF( allocatestatus /= 0 ) stop "*** Not enough memory ***"
149 ALLOCATE ( c( maxm, maxn ), stat = allocatestatus )
150 IF( allocatestatus /= 0 ) stop "*** Not enough memory ***"
151 ALLOCATE ( cc( maxm, maxn ), stat = allocatestatus )
152 IF( allocatestatus /= 0 ) stop "*** Not enough memory ***"
153 ALLOCATE ( x( maxm, maxn ), stat = allocatestatus )
154 IF( allocatestatus /= 0 ) stop "*** Not enough memory ***"
155 ALLOCATE ( swork( ldswork, 103 ), stat = allocatestatus )
156 IF( allocatestatus /= 0 ) stop "*** Not enough memory ***"
157* ..
158* .. Executable Statements ..
159*
160* Get machine parameters
161*
162 eps = dlamch( 'P' )
163 smlnum = dlamch( 'S' ) / eps
164 bignum = one / smlnum
165*
166* Expect INFO = 0
167 vm( 1 ) = one
168* Expect INFO = 1
169 vm( 2 ) = 0.05d+0
170*
171* Begin test loop
172*
173 ninfo( 1 ) = 0
174 ninfo( 2 ) = 0
175 nfail( 1 ) = 0
176 nfail( 2 ) = 0
177 nfail( 3 ) = 0
178 rmax( 1 ) = zero
179 rmax( 2 ) = zero
180 knt = 0
181 iseed( 1 ) = 1
182 iseed( 2 ) = 1
183 iseed( 3 ) = 1
184 iseed( 4 ) = 1
185 scale = one
186 scale3 = one
187 DO j = 1, 2
188 DO isgn = -1, 1, 2
189* Reset seed (overwritten by LATMR)
190 iseed( 1 ) = 1
191 iseed( 2 ) = 1
192 iseed( 3 ) = 1
193 iseed( 4 ) = 1
194 DO m = 32, maxm, 51
195 kla = 0
196 kua = m - 1
197 CALL zlatmr( m, m, 'S', iseed, 'N', d,
198 $ 6, one, cone, 'T', 'N',
199 $ duml, 1, one, dumr, 1, one,
200 $ 'N', iwork, kla, kua, zero,
201 $ one, 'NO', a, maxm, iwork,
202 $ iinfo )
203 DO i = 1, m
204 a( i, i ) = a( i, i ) * vm( j )
205 END DO
206 anrm = zlange( 'M', m, m, a, maxm, dum )
207 DO n = 51, maxn, 47
208 klb = 0
209 kub = n - 1
210 CALL zlatmr( n, n, 'S', iseed, 'N', d,
211 $ 6, one, cone, 'T', 'N',
212 $ duml, 1, one, dumr, 1, one,
213 $ 'N', iwork, klb, kub, zero,
214 $ one, 'NO', b, maxn, iwork,
215 $ iinfo )
216 DO i = 1, n
217 b( i, i ) = b( i, i ) * vm( j )
218 END DO
219 bnrm = zlange( 'M', n, n, b, maxn, dum )
220 tnrm = max( anrm, bnrm )
221 CALL zlatmr( m, n, 'S', iseed, 'N', d,
222 $ 6, one, cone, 'T', 'N',
223 $ duml, 1, one, dumr, 1, one,
224 $ 'N', iwork, m, n, zero, one,
225 $ 'NO', c, maxm, iwork, iinfo )
226 DO itrana = 1, 2
227 IF( itrana.EQ.1 )
228 $ trana = 'N'
229 IF( itrana.EQ.2 )
230 $ trana = 'C'
231 DO itranb = 1, 2
232 IF( itranb.EQ.1 )
233 $ tranb = 'N'
234 IF( itranb.EQ.2 )
235 $ tranb = 'C'
236 knt = knt + 1
237*
238 CALL zlacpy( 'All', m, n, c, maxm, x, maxm)
239 CALL zlacpy( 'All', m, n, c, maxm, cc, maxm)
240 CALL ztrsyl( trana, tranb, isgn, m, n,
241 $ a, maxm, b, maxn, x, maxm,
242 $ scale, iinfo )
243 IF( iinfo.NE.0 )
244 $ ninfo( 1 ) = ninfo( 1 ) + 1
245 xnrm = zlange( 'M', m, n, x, maxm, dum )
246 rmul = cone
247 IF( xnrm.GT.one .AND. tnrm.GT.one ) THEN
248 IF( xnrm.GT.bignum / tnrm ) THEN
249 rmul = cone / max( xnrm, tnrm )
250 END IF
251 END IF
252 CALL zgemm( trana, 'N', m, n, m, rmul,
253 $ a, maxm, x, maxm, -scale*rmul,
254 $ cc, maxm )
255 CALL zgemm( 'N', tranb, m, n, n,
256 $ dble( isgn )*rmul, x, maxm, b,
257 $ maxn, cone, cc, maxm )
258 res1 = zlange( 'M', m, n, cc, maxm, dum )
259 res = res1 / max( smlnum, smlnum*xnrm,
260 $ ( ( abs( rmul )*tnrm )*eps )*xnrm )
261 IF( res.GT.thresh )
262 $ nfail( 1 ) = nfail( 1 ) + 1
263 IF( res.GT.rmax( 1 ) )
264 $ rmax( 1 ) = res
265*
266 CALL zlacpy( 'All', m, n, c, maxm, x, maxm )
267 CALL zlacpy( 'All', m, n, c, maxm, cc, maxm )
268 CALL ztrsyl3( trana, tranb, isgn, m, n,
269 $ a, maxm, b, maxn, x, maxm,
270 $ scale3, swork, ldswork, info)
271 IF( info.NE.0 )
272 $ ninfo( 2 ) = ninfo( 2 ) + 1
273 xnrm = zlange( 'M', m, n, x, maxm, dum )
274 rmul = cone
275 IF( xnrm.GT.one .AND. tnrm.GT.one ) THEN
276 IF( xnrm.GT.bignum / tnrm ) THEN
277 rmul = cone / max( xnrm, tnrm )
278 END IF
279 END IF
280 CALL zgemm( trana, 'N', m, n, m, rmul,
281 $ a, maxm, x, maxm, -scale3*rmul,
282 $ cc, maxm )
283 CALL zgemm( 'N', tranb, m, n, n,
284 $ dble( isgn )*rmul, x, maxm, b,
285 $ maxn, cone, cc, maxm )
286 res1 = zlange( 'M', m, n, cc, maxm, dum )
287 res = res1 / max( smlnum, smlnum*xnrm,
288 $ ( ( abs( rmul )*tnrm )*eps )*xnrm )
289* Verify that TRSYL3 only flushes if TRSYL flushes (but
290* there may be cases where TRSYL3 avoid flushing).
291 IF( scale3.EQ.zero .AND. scale.GT.zero .OR.
292 $ iinfo.NE.info ) THEN
293 nfail( 3 ) = nfail( 3 ) + 1
294 END IF
295 IF( res.GT.thresh .OR. disnan( res ) )
296 $ nfail( 2 ) = nfail( 2 ) + 1
297 IF( res.GT.rmax( 2 ) )
298 $ rmax( 2 ) = res
299 END DO
300 END DO
301 END DO
302 END DO
303 END DO
304 END DO
305*
306 DEALLOCATE (a, stat = allocatestatus)
307 DEALLOCATE (b, stat = allocatestatus)
308 DEALLOCATE (c, stat = allocatestatus)
309 DEALLOCATE (cc, stat = allocatestatus)
310 DEALLOCATE (x, stat = allocatestatus)
311 DEALLOCATE (swork, stat = allocatestatus)
312*
313 RETURN
314*
315* End of ZSYL01
316*
zgemm
subroutine zgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
ZGEMM
Definition zgemm.f:188
disnan
logical function disnan(din)
DISNAN tests input for NaN.
Definition disnan.f:59
zlacpy
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:103
dlamch
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
zlange
double precision function zlange(norm, m, n, a, lda, work)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition zlange.f:115
ztrsyl3
subroutine ztrsyl3(trana, tranb, isgn, m, n, a, lda, b, ldb, c, ldc, scale, swork, ldswork, info)
ZTRSYL3
Definition ztrsyl3.f:157
ztrsyl
subroutine ztrsyl(trana, tranb, isgn, m, n, a, lda, b, ldb, c, ldc, scale, info)
ZTRSYL
Definition ztrsyl.f:157
zlatmr
subroutine zlatmr(m, n, dist, iseed, sym, d, mode, cond, dmax, rsign, grade, dl, model, condl, dr, moder, condr, pivtng, ipivot, kl, ku, sparse, anorm, pack, a, lda, iwork, info)
ZLATMR
Definition zlatmr.f:490
Here is the call graph for this function:
Here is the caller graph for this function: