claqr0.f.html
来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 626 行 · 第 1/4 页
HTML
626 行
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Nibble crossover point ====
</span><span class="comment">*</span><span class="comment">
</span> NIBBLE = <a name="ILAENV.299"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 14, <span class="string">'<a name="CLAQR0.299"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a>'</span>, JBCMPZ, N, ILO, IHI, LWORK )
NIBBLE = MAX( 0, NIBBLE )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Accumulate reflections during ttswp? Use block
</span><span class="comment">*</span><span class="comment"> . 2-by-2 structure during matrix-matrix multiply? ====
</span><span class="comment">*</span><span class="comment">
</span> KACC22 = <a name="ILAENV.305"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 16, <span class="string">'<a name="CLAQR0.305"></a><a href="claqr0.f.html#CLAQR0.1">CLAQR0</a>'</span>, JBCMPZ, N, ILO, IHI, LWORK )
KACC22 = MAX( 0, KACC22 )
KACC22 = MIN( 2, KACC22 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== NWMAX = the largest possible deflation window for
</span><span class="comment">*</span><span class="comment"> . which there is sufficient workspace. ====
</span><span class="comment">*</span><span class="comment">
</span> NWMAX = MIN( ( N-1 ) / 3, LWORK / 2 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== NSMAX = the Largest number of simultaneous shifts
</span><span class="comment">*</span><span class="comment"> . for which there is sufficient workspace. ====
</span><span class="comment">*</span><span class="comment">
</span> NSMAX = MIN( ( N+6 ) / 9, 2*LWORK / 3 )
NSMAX = NSMAX - MOD( NSMAX, 2 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== NDFL: an iteration count restarted at deflation. ====
</span><span class="comment">*</span><span class="comment">
</span> NDFL = 1
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== ITMAX = iteration limit ====
</span><span class="comment">*</span><span class="comment">
</span> ITMAX = MAX( 30, 2*KEXSH )*MAX( 10, ( IHI-ILO+1 ) )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Last row and column in the active block ====
</span><span class="comment">*</span><span class="comment">
</span> KBOT = IHI
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Main Loop ====
</span><span class="comment">*</span><span class="comment">
</span> DO 70 IT = 1, ITMAX
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Done when KBOT falls below ILO ====
</span><span class="comment">*</span><span class="comment">
</span> IF( KBOT.LT.ILO )
$ GO TO 80
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Locate active block ====
</span><span class="comment">*</span><span class="comment">
</span> DO 10 K = KBOT, ILO + 1, -1
IF( H( K, K-1 ).EQ.ZERO )
$ GO TO 20
10 CONTINUE
K = ILO
20 CONTINUE
KTOP = K
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Select deflation window size ====
</span><span class="comment">*</span><span class="comment">
</span> NH = KBOT - KTOP + 1
IF( NDFL.LT.KEXNW .OR. NH.LT.NW ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Typical deflation window. If possible and
</span><span class="comment">*</span><span class="comment"> . advisable, nibble the entire active block.
</span><span class="comment">*</span><span class="comment"> . If not, use size NWR or NWR+1 depending upon
</span><span class="comment">*</span><span class="comment"> . which has the smaller corresponding subdiagonal
</span><span class="comment">*</span><span class="comment"> . entry (a heuristic). ====
</span><span class="comment">*</span><span class="comment">
</span> NWINC = .TRUE.
IF( NH.LE.MIN( NMIN, NWMAX ) ) THEN
NW = NH
ELSE
NW = MIN( NWR, NH, NWMAX )
IF( NW.LT.NWMAX ) THEN
IF( NW.GE.NH-1 ) THEN
NW = NH
ELSE
KWTOP = KBOT - NW + 1
IF( CABS1( H( KWTOP, KWTOP-1 ) ).GT.
$ CABS1( H( KWTOP-1, KWTOP-2 ) ) )NW = NW + 1
END IF
END IF
END IF
ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Exceptional deflation window. If there have
</span><span class="comment">*</span><span class="comment"> . been no deflations in KEXNW or more iterations,
</span><span class="comment">*</span><span class="comment"> . then vary the deflation window size. At first,
</span><span class="comment">*</span><span class="comment"> . because, larger windows are, in general, more
</span><span class="comment">*</span><span class="comment"> . powerful than smaller ones, rapidly increase the
</span><span class="comment">*</span><span class="comment"> . window up to the maximum reasonable and possible.
</span><span class="comment">*</span><span class="comment"> . Then maybe try a slightly smaller window. ====
</span><span class="comment">*</span><span class="comment">
</span> IF( NWINC .AND. NW.LT.MIN( NWMAX, NH ) ) THEN
NW = MIN( NWMAX, NH, 2*NW )
ELSE
NWINC = .FALSE.
IF( NW.EQ.NH .AND. NH.GT.2 )
$ NW = NH - 1
END IF
END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Aggressive early deflation:
</span><span class="comment">*</span><span class="comment"> . split workspace under the subdiagonal into
</span><span class="comment">*</span><span class="comment"> . - an nw-by-nw work array V in the lower
</span><span class="comment">*</span><span class="comment"> . left-hand-corner,
</span><span class="comment">*</span><span class="comment"> . - an NW-by-at-least-NW-but-more-is-better
</span><span class="comment">*</span><span class="comment"> . (NW-by-NHO) horizontal work array along
</span><span class="comment">*</span><span class="comment"> . the bottom edge,
</span><span class="comment">*</span><span class="comment"> . - an at-least-NW-but-more-is-better (NHV-by-NW)
</span><span class="comment">*</span><span class="comment"> . vertical work array along the left-hand-edge.
</span><span class="comment">*</span><span class="comment"> . ====
</span><span class="comment">*</span><span class="comment">
</span> KV = N - NW + 1
KT = NW + 1
NHO = ( N-NW-1 ) - KT + 1
KWV = NW + 2
NVE = ( N-NW ) - KWV + 1
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Aggressive early deflation ====
</span><span class="comment">*</span><span class="comment">
</span> CALL <a name="CLAQR3.415"></a><a href="claqr3.f.html#CLAQR3.1">CLAQR3</a>( WANTT, WANTZ, N, KTOP, KBOT, NW, H, LDH, ILOZ,
$ IHIZ, Z, LDZ, LS, LD, W, H( KV, 1 ), LDH, NHO,
$ H( KV, KT ), LDH, NVE, H( KWV, 1 ), LDH, WORK,
$ LWORK )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Adjust KBOT accounting for new deflations. ====
</span><span class="comment">*</span><span class="comment">
</span> KBOT = KBOT - LD
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== KS points to the shifts. ====
</span><span class="comment">*</span><span class="comment">
</span> KS = KBOT - LS + 1
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== Skip an expensive QR sweep if there is a (partly
</span><span class="comment">*</span><span class="comment"> . heuristic) reason to expect that many eigenvalues
</span><span class="comment">*</span><span class="comment"> . will deflate without it. Here, the QR sweep is
</span><span class="comment">*</span><span class="comment"> . skipped if many eigenvalues have just been deflated
</span><span class="comment">*</span><span class="comment"> . or if the remaining active block is small.
</span><span class="comment">*</span><span class="comment">
</span> IF( ( LD.EQ.0 ) .OR. ( ( 100*LD.LE.NW*NIBBLE ) .AND. ( KBOT-
$ KTOP+1.GT.MIN( NMIN, NWMAX ) ) ) ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== NS = nominal number of simultaneous shifts.
</span><span class="comment">*</span><span class="comment"> . This may be lowered (slightly) if <a name="CLAQR3.438"></a><a href="claqr3.f.html#CLAQR3.1">CLAQR3</a>
</span><span class="comment">*</span><span class="comment"> . did not provide that many shifts. ====
</span><span class="comment">*</span><span class="comment">
</span> NS = MIN( NSMAX, NSR, MAX( 2, KBOT-KTOP ) )
NS = NS - MOD( NS, 2 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> ==== If there have been no deflations
</span><span class="comment">*</span><span class="comment"> . in a multiple of KEXSH iterations,
</span><span class="comment">*</span><span class="comment"> . then try exceptional shifts.
</span><span class="comment">*</span><span class="comment"> . Otherwise use shifts provided by
</span><span class="comment">*</span><span class="comment"> . <a name="CLAQR3.448"></a><a href="claqr3.f.html#CLAQR3.1">CLAQR3</a> above or from the eigenvalues
</span><span class="comment">*</span><span class="comment"> . of a trailing principal submatrix. ====
</span><span class="comment">*</span><span class="comment">
</span> IF( MOD( NDFL, KEXSH ).EQ.0 ) THEN
KS = KBOT - NS + 1
⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?