Annotation of rpl/lapack/lapack/dlaebz.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE DLAEBZ( IJOB, NITMAX, N, MMAX, MINP, NBMIN, ABSTOL,
! 2: $ RELTOL, PIVMIN, D, E, E2, NVAL, AB, C, MOUT,
! 3: $ NAB, WORK, IWORK, INFO )
! 4: *
! 5: * -- LAPACK auxiliary routine (version 3.2) --
! 6: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 7: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 8: * November 2006
! 9: *
! 10: * .. Scalar Arguments ..
! 11: INTEGER IJOB, INFO, MINP, MMAX, MOUT, N, NBMIN, NITMAX
! 12: DOUBLE PRECISION ABSTOL, PIVMIN, RELTOL
! 13: * ..
! 14: * .. Array Arguments ..
! 15: INTEGER IWORK( * ), NAB( MMAX, * ), NVAL( * )
! 16: DOUBLE PRECISION AB( MMAX, * ), C( * ), D( * ), E( * ), E2( * ),
! 17: $ WORK( * )
! 18: * ..
! 19: *
! 20: * Purpose
! 21: * =======
! 22: *
! 23: * DLAEBZ contains the iteration loops which compute and use the
! 24: * function N(w), which is the count of eigenvalues of a symmetric
! 25: * tridiagonal matrix T less than or equal to its argument w. It
! 26: * performs a choice of two types of loops:
! 27: *
! 28: * IJOB=1, followed by
! 29: * IJOB=2: It takes as input a list of intervals and returns a list of
! 30: * sufficiently small intervals whose union contains the same
! 31: * eigenvalues as the union of the original intervals.
! 32: * The input intervals are (AB(j,1),AB(j,2)], j=1,...,MINP.
! 33: * The output interval (AB(j,1),AB(j,2)] will contain
! 34: * eigenvalues NAB(j,1)+1,...,NAB(j,2), where 1 <= j <= MOUT.
! 35: *
! 36: * IJOB=3: It performs a binary search in each input interval
! 37: * (AB(j,1),AB(j,2)] for a point w(j) such that
! 38: * N(w(j))=NVAL(j), and uses C(j) as the starting point of
! 39: * the search. If such a w(j) is found, then on output
! 40: * AB(j,1)=AB(j,2)=w. If no such w(j) is found, then on output
! 41: * (AB(j,1),AB(j,2)] will be a small interval containing the
! 42: * point where N(w) jumps through NVAL(j), unless that point
! 43: * lies outside the initial interval.
! 44: *
! 45: * Note that the intervals are in all cases half-open intervals,
! 46: * i.e., of the form (a,b] , which includes b but not a .
! 47: *
! 48: * To avoid underflow, the matrix should be scaled so that its largest
! 49: * element is no greater than overflow**(1/2) * underflow**(1/4)
! 50: * in absolute value. To assure the most accurate computation
! 51: * of small eigenvalues, the matrix should be scaled to be
! 52: * not much smaller than that, either.
! 53: *
! 54: * See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal
! 55: * Matrix", Report CS41, Computer Science Dept., Stanford
! 56: * University, July 21, 1966
! 57: *
! 58: * Note: the arguments are, in general, *not* checked for unreasonable
! 59: * values.
! 60: *
! 61: * Arguments
! 62: * =========
! 63: *
! 64: * IJOB (input) INTEGER
! 65: * Specifies what is to be done:
! 66: * = 1: Compute NAB for the initial intervals.
! 67: * = 2: Perform bisection iteration to find eigenvalues of T.
! 68: * = 3: Perform bisection iteration to invert N(w), i.e.,
! 69: * to find a point which has a specified number of
! 70: * eigenvalues of T to its left.
! 71: * Other values will cause DLAEBZ to return with INFO=-1.
! 72: *
! 73: * NITMAX (input) INTEGER
! 74: * The maximum number of "levels" of bisection to be
! 75: * performed, i.e., an interval of width W will not be made
! 76: * smaller than 2^(-NITMAX) * W. If not all intervals
! 77: * have converged after NITMAX iterations, then INFO is set
! 78: * to the number of non-converged intervals.
! 79: *
! 80: * N (input) INTEGER
! 81: * The dimension n of the tridiagonal matrix T. It must be at
! 82: * least 1.
! 83: *
! 84: * MMAX (input) INTEGER
! 85: * The maximum number of intervals. If more than MMAX intervals
! 86: * are generated, then DLAEBZ will quit with INFO=MMAX+1.
! 87: *
! 88: * MINP (input) INTEGER
! 89: * The initial number of intervals. It may not be greater than
! 90: * MMAX.
! 91: *
! 92: * NBMIN (input) INTEGER
! 93: * The smallest number of intervals that should be processed
! 94: * using a vector loop. If zero, then only the scalar loop
! 95: * will be used.
! 96: *
! 97: * ABSTOL (input) DOUBLE PRECISION
! 98: * The minimum (absolute) width of an interval. When an
! 99: * interval is narrower than ABSTOL, or than RELTOL times the
! 100: * larger (in magnitude) endpoint, then it is considered to be
! 101: * sufficiently small, i.e., converged. This must be at least
! 102: * zero.
! 103: *
! 104: * RELTOL (input) DOUBLE PRECISION
! 105: * The minimum relative width of an interval. When an interval
! 106: * is narrower than ABSTOL, or than RELTOL times the larger (in
! 107: * magnitude) endpoint, then it is considered to be
! 108: * sufficiently small, i.e., converged. Note: this should
! 109: * always be at least radix*machine epsilon.
! 110: *
! 111: * PIVMIN (input) DOUBLE PRECISION
! 112: * The minimum absolute value of a "pivot" in the Sturm
! 113: * sequence loop. This *must* be at least max |e(j)**2| *
! 114: * safe_min and at least safe_min, where safe_min is at least
! 115: * the smallest number that can divide one without overflow.
! 116: *
! 117: * D (input) DOUBLE PRECISION array, dimension (N)
! 118: * The diagonal elements of the tridiagonal matrix T.
! 119: *
! 120: * E (input) DOUBLE PRECISION array, dimension (N)
! 121: * The offdiagonal elements of the tridiagonal matrix T in
! 122: * positions 1 through N-1. E(N) is arbitrary.
! 123: *
! 124: * E2 (input) DOUBLE PRECISION array, dimension (N)
! 125: * The squares of the offdiagonal elements of the tridiagonal
! 126: * matrix T. E2(N) is ignored.
! 127: *
! 128: * NVAL (input/output) INTEGER array, dimension (MINP)
! 129: * If IJOB=1 or 2, not referenced.
! 130: * If IJOB=3, the desired values of N(w). The elements of NVAL
! 131: * will be reordered to correspond with the intervals in AB.
! 132: * Thus, NVAL(j) on output will not, in general be the same as
! 133: * NVAL(j) on input, but it will correspond with the interval
! 134: * (AB(j,1),AB(j,2)] on output.
! 135: *
! 136: * AB (input/output) DOUBLE PRECISION array, dimension (MMAX,2)
! 137: * The endpoints of the intervals. AB(j,1) is a(j), the left
! 138: * endpoint of the j-th interval, and AB(j,2) is b(j), the
! 139: * right endpoint of the j-th interval. The input intervals
! 140: * will, in general, be modified, split, and reordered by the
! 141: * calculation.
! 142: *
! 143: * C (input/output) DOUBLE PRECISION array, dimension (MMAX)
! 144: * If IJOB=1, ignored.
! 145: * If IJOB=2, workspace.
! 146: * If IJOB=3, then on input C(j) should be initialized to the
! 147: * first search point in the binary search.
! 148: *
! 149: * MOUT (output) INTEGER
! 150: * If IJOB=1, the number of eigenvalues in the intervals.
! 151: * If IJOB=2 or 3, the number of intervals output.
! 152: * If IJOB=3, MOUT will equal MINP.
! 153: *
! 154: * NAB (input/output) INTEGER array, dimension (MMAX,2)
! 155: * If IJOB=1, then on output NAB(i,j) will be set to N(AB(i,j)).
! 156: * If IJOB=2, then on input, NAB(i,j) should be set. It must
! 157: * satisfy the condition:
! 158: * N(AB(i,1)) <= NAB(i,1) <= NAB(i,2) <= N(AB(i,2)),
! 159: * which means that in interval i only eigenvalues
! 160: * NAB(i,1)+1,...,NAB(i,2) will be considered. Usually,
! 161: * NAB(i,j)=N(AB(i,j)), from a previous call to DLAEBZ with
! 162: * IJOB=1.
! 163: * On output, NAB(i,j) will contain
! 164: * max(na(k),min(nb(k),N(AB(i,j)))), where k is the index of
! 165: * the input interval that the output interval
! 166: * (AB(j,1),AB(j,2)] came from, and na(k) and nb(k) are the
! 167: * the input values of NAB(k,1) and NAB(k,2).
! 168: * If IJOB=3, then on output, NAB(i,j) contains N(AB(i,j)),
! 169: * unless N(w) > NVAL(i) for all search points w , in which
! 170: * case NAB(i,1) will not be modified, i.e., the output
! 171: * value will be the same as the input value (modulo
! 172: * reorderings -- see NVAL and AB), or unless N(w) < NVAL(i)
! 173: * for all search points w , in which case NAB(i,2) will
! 174: * not be modified. Normally, NAB should be set to some
! 175: * distinctive value(s) before DLAEBZ is called.
! 176: *
! 177: * WORK (workspace) DOUBLE PRECISION array, dimension (MMAX)
! 178: * Workspace.
! 179: *
! 180: * IWORK (workspace) INTEGER array, dimension (MMAX)
! 181: * Workspace.
! 182: *
! 183: * INFO (output) INTEGER
! 184: * = 0: All intervals converged.
! 185: * = 1--MMAX: The last INFO intervals did not converge.
! 186: * = MMAX+1: More than MMAX intervals were generated.
! 187: *
! 188: * Further Details
! 189: * ===============
! 190: *
! 191: * This routine is intended to be called only by other LAPACK
! 192: * routines, thus the interface is less user-friendly. It is intended
! 193: * for two purposes:
! 194: *
! 195: * (a) finding eigenvalues. In this case, DLAEBZ should have one or
! 196: * more initial intervals set up in AB, and DLAEBZ should be called
! 197: * with IJOB=1. This sets up NAB, and also counts the eigenvalues.
! 198: * Intervals with no eigenvalues would usually be thrown out at
! 199: * this point. Also, if not all the eigenvalues in an interval i
! 200: * are desired, NAB(i,1) can be increased or NAB(i,2) decreased.
! 201: * For example, set NAB(i,1)=NAB(i,2)-1 to get the largest
! 202: * eigenvalue. DLAEBZ is then called with IJOB=2 and MMAX
! 203: * no smaller than the value of MOUT returned by the call with
! 204: * IJOB=1. After this (IJOB=2) call, eigenvalues NAB(i,1)+1
! 205: * through NAB(i,2) are approximately AB(i,1) (or AB(i,2)) to the
! 206: * tolerance specified by ABSTOL and RELTOL.
! 207: *
! 208: * (b) finding an interval (a',b'] containing eigenvalues w(f),...,w(l).
! 209: * In this case, start with a Gershgorin interval (a,b). Set up
! 210: * AB to contain 2 search intervals, both initially (a,b). One
! 211: * NVAL element should contain f-1 and the other should contain l
! 212: * , while C should contain a and b, resp. NAB(i,1) should be -1
! 213: * and NAB(i,2) should be N+1, to flag an error if the desired
! 214: * interval does not lie in (a,b). DLAEBZ is then called with
! 215: * IJOB=3. On exit, if w(f-1) < w(f), then one of the intervals --
! 216: * j -- will have AB(j,1)=AB(j,2) and NAB(j,1)=NAB(j,2)=f-1, while
! 217: * if, to the specified tolerance, w(f-k)=...=w(f+r), k > 0 and r
! 218: * >= 0, then the interval will have N(AB(j,1))=NAB(j,1)=f-k and
! 219: * N(AB(j,2))=NAB(j,2)=f+r. The cases w(l) < w(l+1) and
! 220: * w(l-r)=...=w(l+k) are handled similarly.
! 221: *
! 222: * =====================================================================
! 223: *
! 224: * .. Parameters ..
! 225: DOUBLE PRECISION ZERO, TWO, HALF
! 226: PARAMETER ( ZERO = 0.0D0, TWO = 2.0D0,
! 227: $ HALF = 1.0D0 / TWO )
! 228: * ..
! 229: * .. Local Scalars ..
! 230: INTEGER ITMP1, ITMP2, J, JI, JIT, JP, KF, KFNEW, KL,
! 231: $ KLNEW
! 232: DOUBLE PRECISION TMP1, TMP2
! 233: * ..
! 234: * .. Intrinsic Functions ..
! 235: INTRINSIC ABS, MAX, MIN
! 236: * ..
! 237: * .. Executable Statements ..
! 238: *
! 239: * Check for Errors
! 240: *
! 241: INFO = 0
! 242: IF( IJOB.LT.1 .OR. IJOB.GT.3 ) THEN
! 243: INFO = -1
! 244: RETURN
! 245: END IF
! 246: *
! 247: * Initialize NAB
! 248: *
! 249: IF( IJOB.EQ.1 ) THEN
! 250: *
! 251: * Compute the number of eigenvalues in the initial intervals.
! 252: *
! 253: MOUT = 0
! 254: *DIR$ NOVECTOR
! 255: DO 30 JI = 1, MINP
! 256: DO 20 JP = 1, 2
! 257: TMP1 = D( 1 ) - AB( JI, JP )
! 258: IF( ABS( TMP1 ).LT.PIVMIN )
! 259: $ TMP1 = -PIVMIN
! 260: NAB( JI, JP ) = 0
! 261: IF( TMP1.LE.ZERO )
! 262: $ NAB( JI, JP ) = 1
! 263: *
! 264: DO 10 J = 2, N
! 265: TMP1 = D( J ) - E2( J-1 ) / TMP1 - AB( JI, JP )
! 266: IF( ABS( TMP1 ).LT.PIVMIN )
! 267: $ TMP1 = -PIVMIN
! 268: IF( TMP1.LE.ZERO )
! 269: $ NAB( JI, JP ) = NAB( JI, JP ) + 1
! 270: 10 CONTINUE
! 271: 20 CONTINUE
! 272: MOUT = MOUT + NAB( JI, 2 ) - NAB( JI, 1 )
! 273: 30 CONTINUE
! 274: RETURN
! 275: END IF
! 276: *
! 277: * Initialize for loop
! 278: *
! 279: * KF and KL have the following meaning:
! 280: * Intervals 1,...,KF-1 have converged.
! 281: * Intervals KF,...,KL still need to be refined.
! 282: *
! 283: KF = 1
! 284: KL = MINP
! 285: *
! 286: * If IJOB=2, initialize C.
! 287: * If IJOB=3, use the user-supplied starting point.
! 288: *
! 289: IF( IJOB.EQ.2 ) THEN
! 290: DO 40 JI = 1, MINP
! 291: C( JI ) = HALF*( AB( JI, 1 )+AB( JI, 2 ) )
! 292: 40 CONTINUE
! 293: END IF
! 294: *
! 295: * Iteration loop
! 296: *
! 297: DO 130 JIT = 1, NITMAX
! 298: *
! 299: * Loop over intervals
! 300: *
! 301: IF( KL-KF+1.GE.NBMIN .AND. NBMIN.GT.0 ) THEN
! 302: *
! 303: * Begin of Parallel Version of the loop
! 304: *
! 305: DO 60 JI = KF, KL
! 306: *
! 307: * Compute N(c), the number of eigenvalues less than c
! 308: *
! 309: WORK( JI ) = D( 1 ) - C( JI )
! 310: IWORK( JI ) = 0
! 311: IF( WORK( JI ).LE.PIVMIN ) THEN
! 312: IWORK( JI ) = 1
! 313: WORK( JI ) = MIN( WORK( JI ), -PIVMIN )
! 314: END IF
! 315: *
! 316: DO 50 J = 2, N
! 317: WORK( JI ) = D( J ) - E2( J-1 ) / WORK( JI ) - C( JI )
! 318: IF( WORK( JI ).LE.PIVMIN ) THEN
! 319: IWORK( JI ) = IWORK( JI ) + 1
! 320: WORK( JI ) = MIN( WORK( JI ), -PIVMIN )
! 321: END IF
! 322: 50 CONTINUE
! 323: 60 CONTINUE
! 324: *
! 325: IF( IJOB.LE.2 ) THEN
! 326: *
! 327: * IJOB=2: Choose all intervals containing eigenvalues.
! 328: *
! 329: KLNEW = KL
! 330: DO 70 JI = KF, KL
! 331: *
! 332: * Insure that N(w) is monotone
! 333: *
! 334: IWORK( JI ) = MIN( NAB( JI, 2 ),
! 335: $ MAX( NAB( JI, 1 ), IWORK( JI ) ) )
! 336: *
! 337: * Update the Queue -- add intervals if both halves
! 338: * contain eigenvalues.
! 339: *
! 340: IF( IWORK( JI ).EQ.NAB( JI, 2 ) ) THEN
! 341: *
! 342: * No eigenvalue in the upper interval:
! 343: * just use the lower interval.
! 344: *
! 345: AB( JI, 2 ) = C( JI )
! 346: *
! 347: ELSE IF( IWORK( JI ).EQ.NAB( JI, 1 ) ) THEN
! 348: *
! 349: * No eigenvalue in the lower interval:
! 350: * just use the upper interval.
! 351: *
! 352: AB( JI, 1 ) = C( JI )
! 353: ELSE
! 354: KLNEW = KLNEW + 1
! 355: IF( KLNEW.LE.MMAX ) THEN
! 356: *
! 357: * Eigenvalue in both intervals -- add upper to
! 358: * queue.
! 359: *
! 360: AB( KLNEW, 2 ) = AB( JI, 2 )
! 361: NAB( KLNEW, 2 ) = NAB( JI, 2 )
! 362: AB( KLNEW, 1 ) = C( JI )
! 363: NAB( KLNEW, 1 ) = IWORK( JI )
! 364: AB( JI, 2 ) = C( JI )
! 365: NAB( JI, 2 ) = IWORK( JI )
! 366: ELSE
! 367: INFO = MMAX + 1
! 368: END IF
! 369: END IF
! 370: 70 CONTINUE
! 371: IF( INFO.NE.0 )
! 372: $ RETURN
! 373: KL = KLNEW
! 374: ELSE
! 375: *
! 376: * IJOB=3: Binary search. Keep only the interval containing
! 377: * w s.t. N(w) = NVAL
! 378: *
! 379: DO 80 JI = KF, KL
! 380: IF( IWORK( JI ).LE.NVAL( JI ) ) THEN
! 381: AB( JI, 1 ) = C( JI )
! 382: NAB( JI, 1 ) = IWORK( JI )
! 383: END IF
! 384: IF( IWORK( JI ).GE.NVAL( JI ) ) THEN
! 385: AB( JI, 2 ) = C( JI )
! 386: NAB( JI, 2 ) = IWORK( JI )
! 387: END IF
! 388: 80 CONTINUE
! 389: END IF
! 390: *
! 391: ELSE
! 392: *
! 393: * End of Parallel Version of the loop
! 394: *
! 395: * Begin of Serial Version of the loop
! 396: *
! 397: KLNEW = KL
! 398: DO 100 JI = KF, KL
! 399: *
! 400: * Compute N(w), the number of eigenvalues less than w
! 401: *
! 402: TMP1 = C( JI )
! 403: TMP2 = D( 1 ) - TMP1
! 404: ITMP1 = 0
! 405: IF( TMP2.LE.PIVMIN ) THEN
! 406: ITMP1 = 1
! 407: TMP2 = MIN( TMP2, -PIVMIN )
! 408: END IF
! 409: *
! 410: * A series of compiler directives to defeat vectorization
! 411: * for the next loop
! 412: *
! 413: *$PL$ CMCHAR=' '
! 414: CDIR$ NEXTSCALAR
! 415: C$DIR SCALAR
! 416: CDIR$ NEXT SCALAR
! 417: CVD$L NOVECTOR
! 418: CDEC$ NOVECTOR
! 419: CVD$ NOVECTOR
! 420: *VDIR NOVECTOR
! 421: *VOCL LOOP,SCALAR
! 422: CIBM PREFER SCALAR
! 423: *$PL$ CMCHAR='*'
! 424: *
! 425: DO 90 J = 2, N
! 426: TMP2 = D( J ) - E2( J-1 ) / TMP2 - TMP1
! 427: IF( TMP2.LE.PIVMIN ) THEN
! 428: ITMP1 = ITMP1 + 1
! 429: TMP2 = MIN( TMP2, -PIVMIN )
! 430: END IF
! 431: 90 CONTINUE
! 432: *
! 433: IF( IJOB.LE.2 ) THEN
! 434: *
! 435: * IJOB=2: Choose all intervals containing eigenvalues.
! 436: *
! 437: * Insure that N(w) is monotone
! 438: *
! 439: ITMP1 = MIN( NAB( JI, 2 ),
! 440: $ MAX( NAB( JI, 1 ), ITMP1 ) )
! 441: *
! 442: * Update the Queue -- add intervals if both halves
! 443: * contain eigenvalues.
! 444: *
! 445: IF( ITMP1.EQ.NAB( JI, 2 ) ) THEN
! 446: *
! 447: * No eigenvalue in the upper interval:
! 448: * just use the lower interval.
! 449: *
! 450: AB( JI, 2 ) = TMP1
! 451: *
! 452: ELSE IF( ITMP1.EQ.NAB( JI, 1 ) ) THEN
! 453: *
! 454: * No eigenvalue in the lower interval:
! 455: * just use the upper interval.
! 456: *
! 457: AB( JI, 1 ) = TMP1
! 458: ELSE IF( KLNEW.LT.MMAX ) THEN
! 459: *
! 460: * Eigenvalue in both intervals -- add upper to queue.
! 461: *
! 462: KLNEW = KLNEW + 1
! 463: AB( KLNEW, 2 ) = AB( JI, 2 )
! 464: NAB( KLNEW, 2 ) = NAB( JI, 2 )
! 465: AB( KLNEW, 1 ) = TMP1
! 466: NAB( KLNEW, 1 ) = ITMP1
! 467: AB( JI, 2 ) = TMP1
! 468: NAB( JI, 2 ) = ITMP1
! 469: ELSE
! 470: INFO = MMAX + 1
! 471: RETURN
! 472: END IF
! 473: ELSE
! 474: *
! 475: * IJOB=3: Binary search. Keep only the interval
! 476: * containing w s.t. N(w) = NVAL
! 477: *
! 478: IF( ITMP1.LE.NVAL( JI ) ) THEN
! 479: AB( JI, 1 ) = TMP1
! 480: NAB( JI, 1 ) = ITMP1
! 481: END IF
! 482: IF( ITMP1.GE.NVAL( JI ) ) THEN
! 483: AB( JI, 2 ) = TMP1
! 484: NAB( JI, 2 ) = ITMP1
! 485: END IF
! 486: END IF
! 487: 100 CONTINUE
! 488: KL = KLNEW
! 489: *
! 490: * End of Serial Version of the loop
! 491: *
! 492: END IF
! 493: *
! 494: * Check for convergence
! 495: *
! 496: KFNEW = KF
! 497: DO 110 JI = KF, KL
! 498: TMP1 = ABS( AB( JI, 2 )-AB( JI, 1 ) )
! 499: TMP2 = MAX( ABS( AB( JI, 2 ) ), ABS( AB( JI, 1 ) ) )
! 500: IF( TMP1.LT.MAX( ABSTOL, PIVMIN, RELTOL*TMP2 ) .OR.
! 501: $ NAB( JI, 1 ).GE.NAB( JI, 2 ) ) THEN
! 502: *
! 503: * Converged -- Swap with position KFNEW,
! 504: * then increment KFNEW
! 505: *
! 506: IF( JI.GT.KFNEW ) THEN
! 507: TMP1 = AB( JI, 1 )
! 508: TMP2 = AB( JI, 2 )
! 509: ITMP1 = NAB( JI, 1 )
! 510: ITMP2 = NAB( JI, 2 )
! 511: AB( JI, 1 ) = AB( KFNEW, 1 )
! 512: AB( JI, 2 ) = AB( KFNEW, 2 )
! 513: NAB( JI, 1 ) = NAB( KFNEW, 1 )
! 514: NAB( JI, 2 ) = NAB( KFNEW, 2 )
! 515: AB( KFNEW, 1 ) = TMP1
! 516: AB( KFNEW, 2 ) = TMP2
! 517: NAB( KFNEW, 1 ) = ITMP1
! 518: NAB( KFNEW, 2 ) = ITMP2
! 519: IF( IJOB.EQ.3 ) THEN
! 520: ITMP1 = NVAL( JI )
! 521: NVAL( JI ) = NVAL( KFNEW )
! 522: NVAL( KFNEW ) = ITMP1
! 523: END IF
! 524: END IF
! 525: KFNEW = KFNEW + 1
! 526: END IF
! 527: 110 CONTINUE
! 528: KF = KFNEW
! 529: *
! 530: * Choose Midpoints
! 531: *
! 532: DO 120 JI = KF, KL
! 533: C( JI ) = HALF*( AB( JI, 1 )+AB( JI, 2 ) )
! 534: 120 CONTINUE
! 535: *
! 536: * If no more intervals to refine, quit.
! 537: *
! 538: IF( KF.GT.KL )
! 539: $ GO TO 140
! 540: 130 CONTINUE
! 541: *
! 542: * Converged
! 543: *
! 544: 140 CONTINUE
! 545: INFO = MAX( KL+1-KF, 0 )
! 546: MOUT = KL
! 547: *
! 548: RETURN
! 549: *
! 550: * End of DLAEBZ
! 551: *
! 552: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dlaebz.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack