slarrb.f.html

来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 321 行 · 第 1/2 页

HTML
321
字号
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
 <head>
  <title>slarrb.f</title>
 <meta name="generator" content="emacs 21.3.1; htmlfontify 0.20">
<style type="text/css"><!-- 
body { background: rgb(255, 255, 255);  color: rgb(0, 0, 0);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: none; }
span.default   { background: rgb(255, 255, 255);  color: rgb(0, 0, 0);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: none; }
span.default a { background: rgb(255, 255, 255);  color: rgb(0, 0, 0);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: underline; }
span.comment   { color: rgb(178, 34, 34);  background: rgb(255, 255, 255);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: none; }
span.comment a { color: rgb(178, 34, 34);  background: rgb(255, 255, 255);  font-style: normal;  font-weight: 500;  font-stretch: normal;  font-family: adobe-courier;  font-size: 11pt;  text-decoration: underline; }
 --></style>

 </head>
  <body>

<pre>
      SUBROUTINE <a name="SLARRB.1"></a><a href="slarrb.f.html#SLARRB.1">SLARRB</a>( N, D, LLD, IFIRST, ILAST, RTOL1,
     $                   RTOL2, OFFSET, W, WGAP, WERR, WORK, IWORK,
     $                   PIVMIN, SPDIAM, TWIST, INFO )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  -- LAPACK auxiliary routine (version 3.1) --
</span><span class="comment">*</span><span class="comment">     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
</span><span class="comment">*</span><span class="comment">     November 2006
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     .. Scalar Arguments ..
</span>      INTEGER            IFIRST, ILAST, INFO, N, OFFSET, TWIST
      REAL               PIVMIN, RTOL1, RTOL2, SPDIAM
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Array Arguments ..
</span>      INTEGER            IWORK( * )
      REAL               D( * ), LLD( * ), W( * ),
     $                   WERR( * ), WGAP( * ), WORK( * )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Purpose
</span><span class="comment">*</span><span class="comment">  =======
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Given the relatively robust representation(RRR) L D L^T, <a name="SLARRB.22"></a><a href="slarrb.f.html#SLARRB.1">SLARRB</a>
</span><span class="comment">*</span><span class="comment">  does &quot;limited&quot; bisection to refine the eigenvalues of L D L^T,
</span><span class="comment">*</span><span class="comment">  W( IFIRST-OFFSET ) through W( ILAST-OFFSET ), to more accuracy. Initial
</span><span class="comment">*</span><span class="comment">  guesses for these eigenvalues are input in W, the corresponding estimate
</span><span class="comment">*</span><span class="comment">  of the error in these guesses and their gaps are input in WERR
</span><span class="comment">*</span><span class="comment">  and WGAP, respectively. During bisection, intervals
</span><span class="comment">*</span><span class="comment">  [left, right] are maintained by storing their mid-points and
</span><span class="comment">*</span><span class="comment">  semi-widths in the arrays W and WERR respectively.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Arguments
</span><span class="comment">*</span><span class="comment">  =========
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  N       (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The order of the matrix.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  D       (input) REAL             array, dimension (N)
</span><span class="comment">*</span><span class="comment">          The N diagonal elements of the diagonal matrix D.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  LLD     (input) REAL             array, dimension (N-1)
</span><span class="comment">*</span><span class="comment">          The (N-1) elements L(i)*L(i)*D(i).
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  IFIRST  (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The index of the first eigenvalue to be computed.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  ILAST   (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The index of the last eigenvalue to be computed.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  RTOL1   (input) REAL            
</span><span class="comment">*</span><span class="comment">  RTOL2   (input) REAL            
</span><span class="comment">*</span><span class="comment">          Tolerance for the convergence of the bisection intervals.
</span><span class="comment">*</span><span class="comment">          An interval [LEFT,RIGHT] has converged if
</span><span class="comment">*</span><span class="comment">          RIGHT-LEFT.LT.MAX( RTOL1*GAP, RTOL2*MAX(|LEFT|,|RIGHT|) )
</span><span class="comment">*</span><span class="comment">          where GAP is the (estimated) distance to the nearest
</span><span class="comment">*</span><span class="comment">          eigenvalue.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  OFFSET  (input) INTEGER
</span><span class="comment">*</span><span class="comment">          Offset for the arrays W, WGAP and WERR, i.e., the IFIRST-OFFSET
</span><span class="comment">*</span><span class="comment">          through ILAST-OFFSET elements of these arrays are to be used.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  W       (input/output) REAL             array, dimension (N)
</span><span class="comment">*</span><span class="comment">          On input, W( IFIRST-OFFSET ) through W( ILAST-OFFSET ) are
</span><span class="comment">*</span><span class="comment">          estimates of the eigenvalues of L D L^T indexed IFIRST throug
</span><span class="comment">*</span><span class="comment">          ILAST.
</span><span class="comment">*</span><span class="comment">          On output, these estimates are refined.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  WGAP    (input/output) REAL             array, dimension (N-1)
</span><span class="comment">*</span><span class="comment">          On input, the (estimated) gaps between consecutive
</span><span class="comment">*</span><span class="comment">          eigenvalues of L D L^T, i.e., WGAP(I-OFFSET) is the gap between
</span><span class="comment">*</span><span class="comment">          eigenvalues I and I+1. Note that if IFIRST.EQ.ILAST
</span><span class="comment">*</span><span class="comment">          then WGAP(IFIRST-OFFSET) must be set to ZERO.
</span><span class="comment">*</span><span class="comment">          On output, these gaps are refined.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  WERR    (input/output) REAL             array, dimension (N)
</span><span class="comment">*</span><span class="comment">          On input, WERR( IFIRST-OFFSET ) through WERR( ILAST-OFFSET ) are
</span><span class="comment">*</span><span class="comment">          the errors in the estimates of the corresponding elements in W.
</span><span class="comment">*</span><span class="comment">          On output, these errors are refined.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  WORK    (workspace) REAL             array, dimension (2*N)
</span><span class="comment">*</span><span class="comment">          Workspace.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  IWORK   (workspace) INTEGER array, dimension (2*N)
</span><span class="comment">*</span><span class="comment">          Workspace.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  PIVMIN  (input) DOUBLE PRECISION
</span><span class="comment">*</span><span class="comment">          The minimum pivot in the Sturm sequence.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  SPDIAM  (input) DOUBLE PRECISION
</span><span class="comment">*</span><span class="comment">          The spectral diameter of the matrix.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  TWIST   (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The twist index for the twisted factorization that is used
</span><span class="comment">*</span><span class="comment">          for the negcount.
</span><span class="comment">*</span><span class="comment">          TWIST = N: Compute negcount from L D L^T - LAMBDA I = L+ D+ L+^T
</span><span class="comment">*</span><span class="comment">          TWIST = 1: Compute negcount from L D L^T - LAMBDA I = U- D- U-^T
</span><span class="comment">*</span><span class="comment">          TWIST = R: Compute negcount from L D L^T - LAMBDA I = N(r) D(r) N(r)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  INFO    (output) INTEGER
</span><span class="comment">*</span><span class="comment">          Error flag.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Further Details
</span><span class="comment">*</span><span class="comment">  ===============
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Based on contributions by
</span><span class="comment">*</span><span class="comment">     Beresford Parlett, University of California, Berkeley, USA
</span><span class="comment">*</span><span class="comment">     Jim Demmel, University of California, Berkeley, USA
</span><span class="comment">*</span><span class="comment">     Inderjit Dhillon, University of Texas, Austin, USA
</span><span class="comment">*</span><span class="comment">     Osni Marques, LBNL/NERSC, USA
</span><span class="comment">*</span><span class="comment">     Christof Voemel, University of California, Berkeley, USA
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  =====================================================================
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     .. Parameters ..
</span>      REAL               ZERO, TWO, HALF
      PARAMETER        ( ZERO = 0.0E0, TWO = 2.0E0,
     $                   HALF = 0.5E0 )
      INTEGER   MAXITR
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Scalars ..
</span>      INTEGER            I, I1, II, IP, ITER, K, NEGCNT, NEXT, NINT,
     $                   OLNINT, PREV, R
      REAL               BACK, CVRGD, GAP, LEFT, LGAP, MID, MNWDTH,
     $                   RGAP, RIGHT, TMP, WIDTH
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Functions ..
</span>      INTEGER            <a name="SLANEG.126"></a><a href="slaneg.f.html#SLANEG.1">SLANEG</a>
      EXTERNAL           <a name="SLANEG.127"></a><a href="slaneg.f.html#SLANEG.1">SLANEG</a>
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Intrinsic Functions ..
</span>      INTRINSIC          ABS, MAX, MIN
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Executable Statements ..
</span><span class="comment">*</span><span class="comment">
</span>      INFO = 0
<span class="comment">*</span><span class="comment">
</span>      MAXITR = INT( ( LOG( SPDIAM+PIVMIN )-LOG( PIVMIN ) ) /
     $           LOG( TWO ) ) + 2
      MNWDTH = TWO * PIVMIN
<span class="comment">*</span><span class="comment">
</span>      R = TWIST
      IF((R.LT.1).OR.(R.GT.N)) R = N
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Initialize unconverged intervals in [ WORK(2*I-1), WORK(2*I) ].

⌨️ 快捷键说明

复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?