[BACK]Return to dhgeqz.f CVS log [TXT] [DIR] Up to [local] / rpl / lapack / lapack

Diff for /rpl/lapack/lapack/dhgeqz.f between versions 1.8 and 1.9

version 1.8, 2011/07/22 07:38:06 version 1.9, 2011/11/21 20:42:53
Line 1 Line 1
   *> \brief \b DHGEQZ
   *
   *  =========== DOCUMENTATION ===========
   *
   * Online html documentation available at 
   *            http://www.netlib.org/lapack/explore-html/ 
   *
   *> \htmlonly
   *> Download DHGEQZ + dependencies 
   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dhgeqz.f"> 
   *> [TGZ]</a> 
   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dhgeqz.f"> 
   *> [ZIP]</a> 
   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dhgeqz.f"> 
   *> [TXT]</a>
   *> \endhtmlonly 
   *
   *  Definition:
   *  ===========
   *
   *       SUBROUTINE DHGEQZ( JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT,
   *                          ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK,
   *                          LWORK, INFO )
   * 
   *       .. Scalar Arguments ..
   *       CHARACTER          COMPQ, COMPZ, JOB
   *       INTEGER            IHI, ILO, INFO, LDH, LDQ, LDT, LDZ, LWORK, N
   *       ..
   *       .. Array Arguments ..
   *       DOUBLE PRECISION   ALPHAI( * ), ALPHAR( * ), BETA( * ),
   *      $                   H( LDH, * ), Q( LDQ, * ), T( LDT, * ),
   *      $                   WORK( * ), Z( LDZ, * )
   *       ..
   *  
   *
   *> \par Purpose:
   *  =============
   *>
   *> \verbatim
   *>
   *> DHGEQZ computes the eigenvalues of a real matrix pair (H,T),
   *> where H is an upper Hessenberg matrix and T is upper triangular,
   *> using the double-shift QZ method.
   *> Matrix pairs of this type are produced by the reduction to
   *> generalized upper Hessenberg form of a real matrix pair (A,B):
   *>
   *>    A = Q1*H*Z1**T,  B = Q1*T*Z1**T,
   *>
   *> as computed by DGGHRD.
   *>
   *> If JOB='S', then the Hessenberg-triangular pair (H,T) is
   *> also reduced to generalized Schur form,
   *> 
   *>    H = Q*S*Z**T,  T = Q*P*Z**T,
   *> 
   *> where Q and Z are orthogonal matrices, P is an upper triangular
   *> matrix, and S is a quasi-triangular matrix with 1-by-1 and 2-by-2
   *> diagonal blocks.
   *>
   *> The 1-by-1 blocks correspond to real eigenvalues of the matrix pair
   *> (H,T) and the 2-by-2 blocks correspond to complex conjugate pairs of
   *> eigenvalues.
   *>
   *> Additionally, the 2-by-2 upper triangular diagonal blocks of P
   *> corresponding to 2-by-2 blocks of S are reduced to positive diagonal
   *> form, i.e., if S(j+1,j) is non-zero, then P(j+1,j) = P(j,j+1) = 0,
   *> P(j,j) > 0, and P(j+1,j+1) > 0.
   *>
   *> Optionally, the orthogonal matrix Q from the generalized Schur
   *> factorization may be postmultiplied into an input matrix Q1, and the
   *> orthogonal matrix Z may be postmultiplied into an input matrix Z1.
   *> If Q1 and Z1 are the orthogonal matrices from DGGHRD that reduced
   *> the matrix pair (A,B) to generalized upper Hessenberg form, then the
   *> output matrices Q1*Q and Z1*Z are the orthogonal factors from the
   *> generalized Schur factorization of (A,B):
   *>
   *>    A = (Q1*Q)*S*(Z1*Z)**T,  B = (Q1*Q)*P*(Z1*Z)**T.
   *> 
   *> To avoid overflow, eigenvalues of the matrix pair (H,T) (equivalently,
   *> of (A,B)) are computed as a pair of values (alpha,beta), where alpha is
   *> complex and beta real.
   *> If beta is nonzero, lambda = alpha / beta is an eigenvalue of the
   *> generalized nonsymmetric eigenvalue problem (GNEP)
   *>    A*x = lambda*B*x
   *> and if alpha is nonzero, mu = beta / alpha is an eigenvalue of the
   *> alternate form of the GNEP
   *>    mu*A*y = B*y.
   *> Real eigenvalues can be read directly from the generalized Schur
   *> form: 
   *>   alpha = S(i,i), beta = P(i,i).
   *>
   *> Ref: C.B. Moler & G.W. Stewart, "An Algorithm for Generalized Matrix
   *>      Eigenvalue Problems", SIAM J. Numer. Anal., 10(1973),
   *>      pp. 241--256.
   *> \endverbatim
   *
   *  Arguments:
   *  ==========
   *
   *> \param[in] JOB
   *> \verbatim
   *>          JOB is CHARACTER*1
   *>          = 'E': Compute eigenvalues only;
   *>          = 'S': Compute eigenvalues and the Schur form. 
   *> \endverbatim
   *>
   *> \param[in] COMPQ
   *> \verbatim
   *>          COMPQ is CHARACTER*1
   *>          = 'N': Left Schur vectors (Q) are not computed;
   *>          = 'I': Q is initialized to the unit matrix and the matrix Q
   *>                 of left Schur vectors of (H,T) is returned;
   *>          = 'V': Q must contain an orthogonal matrix Q1 on entry and
   *>                 the product Q1*Q is returned.
   *> \endverbatim
   *>
   *> \param[in] COMPZ
   *> \verbatim
   *>          COMPZ is CHARACTER*1
   *>          = 'N': Right Schur vectors (Z) are not computed;
   *>          = 'I': Z is initialized to the unit matrix and the matrix Z
   *>                 of right Schur vectors of (H,T) is returned;
   *>          = 'V': Z must contain an orthogonal matrix Z1 on entry and
   *>                 the product Z1*Z is returned.
   *> \endverbatim
   *>
   *> \param[in] N
   *> \verbatim
   *>          N is INTEGER
   *>          The order of the matrices H, T, Q, and Z.  N >= 0.
   *> \endverbatim
   *>
   *> \param[in] ILO
   *> \verbatim
   *>          ILO is INTEGER
   *> \endverbatim
   *>
   *> \param[in] IHI
   *> \verbatim
   *>          IHI is INTEGER
   *>          ILO and IHI mark the rows and columns of H which are in
   *>          Hessenberg form.  It is assumed that A is already upper
   *>          triangular in rows and columns 1:ILO-1 and IHI+1:N.
   *>          If N > 0, 1 <= ILO <= IHI <= N; if N = 0, ILO=1 and IHI=0.
   *> \endverbatim
   *>
   *> \param[in,out] H
   *> \verbatim
   *>          H is DOUBLE PRECISION array, dimension (LDH, N)
   *>          On entry, the N-by-N upper Hessenberg matrix H.
   *>          On exit, if JOB = 'S', H contains the upper quasi-triangular
   *>          matrix S from the generalized Schur factorization.
   *>          If JOB = 'E', the diagonal blocks of H match those of S, but
   *>          the rest of H is unspecified.
   *> \endverbatim
   *>
   *> \param[in] LDH
   *> \verbatim
   *>          LDH is INTEGER
   *>          The leading dimension of the array H.  LDH >= max( 1, N ).
   *> \endverbatim
   *>
   *> \param[in,out] T
   *> \verbatim
   *>          T is DOUBLE PRECISION array, dimension (LDT, N)
   *>          On entry, the N-by-N upper triangular matrix T.
   *>          On exit, if JOB = 'S', T contains the upper triangular
   *>          matrix P from the generalized Schur factorization;
   *>          2-by-2 diagonal blocks of P corresponding to 2-by-2 blocks of S
   *>          are reduced to positive diagonal form, i.e., if H(j+1,j) is
   *>          non-zero, then T(j+1,j) = T(j,j+1) = 0, T(j,j) > 0, and
   *>          T(j+1,j+1) > 0.
   *>          If JOB = 'E', the diagonal blocks of T match those of P, but
   *>          the rest of T is unspecified.
   *> \endverbatim
   *>
   *> \param[in] LDT
   *> \verbatim
   *>          LDT is INTEGER
   *>          The leading dimension of the array T.  LDT >= max( 1, N ).
   *> \endverbatim
   *>
   *> \param[out] ALPHAR
   *> \verbatim
   *>          ALPHAR is DOUBLE PRECISION array, dimension (N)
   *>          The real parts of each scalar alpha defining an eigenvalue
   *>          of GNEP.
   *> \endverbatim
   *>
   *> \param[out] ALPHAI
   *> \verbatim
   *>          ALPHAI is DOUBLE PRECISION array, dimension (N)
   *>          The imaginary parts of each scalar alpha defining an
   *>          eigenvalue of GNEP.
   *>          If ALPHAI(j) is zero, then the j-th eigenvalue is real; if
   *>          positive, then the j-th and (j+1)-st eigenvalues are a
   *>          complex conjugate pair, with ALPHAI(j+1) = -ALPHAI(j).
   *> \endverbatim
   *>
   *> \param[out] BETA
   *> \verbatim
   *>          BETA is DOUBLE PRECISION array, dimension (N)
   *>          The scalars beta that define the eigenvalues of GNEP.
   *>          Together, the quantities alpha = (ALPHAR(j),ALPHAI(j)) and
   *>          beta = BETA(j) represent the j-th eigenvalue of the matrix
   *>          pair (A,B), in one of the forms lambda = alpha/beta or
   *>          mu = beta/alpha.  Since either lambda or mu may overflow,
   *>          they should not, in general, be computed.
   *> \endverbatim
   *>
   *> \param[in,out] Q
   *> \verbatim
   *>          Q is DOUBLE PRECISION array, dimension (LDQ, N)
   *>          On entry, if COMPZ = 'V', the orthogonal matrix Q1 used in
   *>          the reduction of (A,B) to generalized Hessenberg form.
   *>          On exit, if COMPZ = 'I', the orthogonal matrix of left Schur
   *>          vectors of (H,T), and if COMPZ = 'V', the orthogonal matrix
   *>          of left Schur vectors of (A,B).
   *>          Not referenced if COMPZ = 'N'.
   *> \endverbatim
   *>
   *> \param[in] LDQ
   *> \verbatim
   *>          LDQ is INTEGER
   *>          The leading dimension of the array Q.  LDQ >= 1.
   *>          If COMPQ='V' or 'I', then LDQ >= N.
   *> \endverbatim
   *>
   *> \param[in,out] Z
   *> \verbatim
   *>          Z is DOUBLE PRECISION array, dimension (LDZ, N)
   *>          On entry, if COMPZ = 'V', the orthogonal matrix Z1 used in
   *>          the reduction of (A,B) to generalized Hessenberg form.
   *>          On exit, if COMPZ = 'I', the orthogonal matrix of
   *>          right Schur vectors of (H,T), and if COMPZ = 'V', the
   *>          orthogonal matrix of right Schur vectors of (A,B).
   *>          Not referenced if COMPZ = 'N'.
   *> \endverbatim
   *>
   *> \param[in] LDZ
   *> \verbatim
   *>          LDZ is INTEGER
   *>          The leading dimension of the array Z.  LDZ >= 1.
   *>          If COMPZ='V' or 'I', then LDZ >= N.
   *> \endverbatim
   *>
   *> \param[out] WORK
   *> \verbatim
   *>          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
   *>          On exit, if INFO >= 0, WORK(1) returns the optimal LWORK.
   *> \endverbatim
   *>
   *> \param[in] LWORK
   *> \verbatim
   *>          LWORK is INTEGER
   *>          The dimension of the array WORK.  LWORK >= max(1,N).
   *>
   *>          If LWORK = -1, then a workspace query is assumed; the routine
   *>          only calculates the optimal size of the WORK array, returns
   *>          this value as the first entry of the WORK array, and no error
   *>          message related to LWORK is issued by XERBLA.
   *> \endverbatim
   *>
   *> \param[out] INFO
   *> \verbatim
   *>          INFO is INTEGER
   *>          = 0: successful exit
   *>          < 0: if INFO = -i, the i-th argument had an illegal value
   *>          = 1,...,N: the QZ iteration did not converge.  (H,T) is not
   *>                     in Schur form, but ALPHAR(i), ALPHAI(i), and
   *>                     BETA(i), i=INFO+1,...,N should be correct.
   *>          = N+1,...,2*N: the shift calculation failed.  (H,T) is not
   *>                     in Schur form, but ALPHAR(i), ALPHAI(i), and
   *>                     BETA(i), i=INFO-N+1,...,N should be correct.
   *> \endverbatim
   *
   *  Authors:
   *  ========
   *
   *> \author Univ. of Tennessee 
   *> \author Univ. of California Berkeley 
   *> \author Univ. of Colorado Denver 
   *> \author NAG Ltd. 
   *
   *> \date November 2011
   *
   *> \ingroup doubleGEcomputational
   *
   *> \par Further Details:
   *  =====================
   *>
   *> \verbatim
   *>
   *>  Iteration counters:
   *>
   *>  JITER  -- counts iterations.
   *>  IITER  -- counts iterations run since ILAST was last
   *>            changed.  This is therefore reset only when a 1-by-1 or
   *>            2-by-2 block deflates off the bottom.
   *> \endverbatim
   *>
   *  =====================================================================
       SUBROUTINE DHGEQZ( JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT,        SUBROUTINE DHGEQZ( JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT,
      $                   ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK,       $                   ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK,
      $                   LWORK, INFO )       $                   LWORK, INFO )
 *  *
 *  -- LAPACK routine (version 3.3.1)                                  --  *  -- LAPACK computational routine (version 3.4.0) --
 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --  *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--  *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 *  -- April 2009                                                      --  *     November 2011
 *  *
 *     .. Scalar Arguments ..  *     .. Scalar Arguments ..
       CHARACTER          COMPQ, COMPZ, JOB        CHARACTER          COMPQ, COMPZ, JOB
Line 17 Line 319
      $                   WORK( * ), Z( LDZ, * )       $                   WORK( * ), Z( LDZ, * )
 *     ..  *     ..
 *  *
 *  Purpose  
 *  =======  
 *  
 *  DHGEQZ computes the eigenvalues of a real matrix pair (H,T),  
 *  where H is an upper Hessenberg matrix and T is upper triangular,  
 *  using the double-shift QZ method.  
 *  Matrix pairs of this type are produced by the reduction to  
 *  generalized upper Hessenberg form of a real matrix pair (A,B):  
 *  
 *     A = Q1*H*Z1**T,  B = Q1*T*Z1**T,  
 *  
 *  as computed by DGGHRD.  
 *  
 *  If JOB='S', then the Hessenberg-triangular pair (H,T) is  
 *  also reduced to generalized Schur form,  
 *    
 *     H = Q*S*Z**T,  T = Q*P*Z**T,  
 *    
 *  where Q and Z are orthogonal matrices, P is an upper triangular  
 *  matrix, and S is a quasi-triangular matrix with 1-by-1 and 2-by-2  
 *  diagonal blocks.  
 *  
 *  The 1-by-1 blocks correspond to real eigenvalues of the matrix pair  
 *  (H,T) and the 2-by-2 blocks correspond to complex conjugate pairs of  
 *  eigenvalues.  
 *  
 *  Additionally, the 2-by-2 upper triangular diagonal blocks of P  
 *  corresponding to 2-by-2 blocks of S are reduced to positive diagonal  
 *  form, i.e., if S(j+1,j) is non-zero, then P(j+1,j) = P(j,j+1) = 0,  
 *  P(j,j) > 0, and P(j+1,j+1) > 0.  
 *  
 *  Optionally, the orthogonal matrix Q from the generalized Schur  
 *  factorization may be postmultiplied into an input matrix Q1, and the  
 *  orthogonal matrix Z may be postmultiplied into an input matrix Z1.  
 *  If Q1 and Z1 are the orthogonal matrices from DGGHRD that reduced  
 *  the matrix pair (A,B) to generalized upper Hessenberg form, then the  
 *  output matrices Q1*Q and Z1*Z are the orthogonal factors from the  
 *  generalized Schur factorization of (A,B):  
 *  
 *     A = (Q1*Q)*S*(Z1*Z)**T,  B = (Q1*Q)*P*(Z1*Z)**T.  
 *    
 *  To avoid overflow, eigenvalues of the matrix pair (H,T) (equivalently,  
 *  of (A,B)) are computed as a pair of values (alpha,beta), where alpha is  
 *  complex and beta real.  
 *  If beta is nonzero, lambda = alpha / beta is an eigenvalue of the  
 *  generalized nonsymmetric eigenvalue problem (GNEP)  
 *     A*x = lambda*B*x  
 *  and if alpha is nonzero, mu = beta / alpha is an eigenvalue of the  
 *  alternate form of the GNEP  
 *     mu*A*y = B*y.  
 *  Real eigenvalues can be read directly from the generalized Schur  
 *  form:   
 *    alpha = S(i,i), beta = P(i,i).  
 *  
 *  Ref: C.B. Moler & G.W. Stewart, "An Algorithm for Generalized Matrix  
 *       Eigenvalue Problems", SIAM J. Numer. Anal., 10(1973),  
 *       pp. 241--256.  
 *  
 *  Arguments  
 *  =========  
 *  
 *  JOB     (input) CHARACTER*1  
 *          = 'E': Compute eigenvalues only;  
 *          = 'S': Compute eigenvalues and the Schur form.   
 *  
 *  COMPQ   (input) CHARACTER*1  
 *          = 'N': Left Schur vectors (Q) are not computed;  
 *          = 'I': Q is initialized to the unit matrix and the matrix Q  
 *                 of left Schur vectors of (H,T) is returned;  
 *          = 'V': Q must contain an orthogonal matrix Q1 on entry and  
 *                 the product Q1*Q is returned.  
 *  
 *  COMPZ   (input) CHARACTER*1  
 *          = 'N': Right Schur vectors (Z) are not computed;  
 *          = 'I': Z is initialized to the unit matrix and the matrix Z  
 *                 of right Schur vectors of (H,T) is returned;  
 *          = 'V': Z must contain an orthogonal matrix Z1 on entry and  
 *                 the product Z1*Z is returned.  
 *  
 *  N       (input) INTEGER  
 *          The order of the matrices H, T, Q, and Z.  N >= 0.  
 *  
 *  ILO     (input) INTEGER  
 *  IHI     (input) INTEGER  
 *          ILO and IHI mark the rows and columns of H which are in  
 *          Hessenberg form.  It is assumed that A is already upper  
 *          triangular in rows and columns 1:ILO-1 and IHI+1:N.  
 *          If N > 0, 1 <= ILO <= IHI <= N; if N = 0, ILO=1 and IHI=0.  
 *  
 *  H       (input/output) DOUBLE PRECISION array, dimension (LDH, N)  
 *          On entry, the N-by-N upper Hessenberg matrix H.  
 *          On exit, if JOB = 'S', H contains the upper quasi-triangular  
 *          matrix S from the generalized Schur factorization.  
 *          If JOB = 'E', the diagonal blocks of H match those of S, but  
 *          the rest of H is unspecified.  
 *  
 *  LDH     (input) INTEGER  
 *          The leading dimension of the array H.  LDH >= max( 1, N ).  
 *  
 *  T       (input/output) DOUBLE PRECISION array, dimension (LDT, N)  
 *          On entry, the N-by-N upper triangular matrix T.  
 *          On exit, if JOB = 'S', T contains the upper triangular  
 *          matrix P from the generalized Schur factorization;  
 *          2-by-2 diagonal blocks of P corresponding to 2-by-2 blocks of S  
 *          are reduced to positive diagonal form, i.e., if H(j+1,j) is  
 *          non-zero, then T(j+1,j) = T(j,j+1) = 0, T(j,j) > 0, and  
 *          T(j+1,j+1) > 0.  
 *          If JOB = 'E', the diagonal blocks of T match those of P, but  
 *          the rest of T is unspecified.  
 *  
 *  LDT     (input) INTEGER  
 *          The leading dimension of the array T.  LDT >= max( 1, N ).  
 *  
 *  ALPHAR  (output) DOUBLE PRECISION array, dimension (N)  
 *          The real parts of each scalar alpha defining an eigenvalue  
 *          of GNEP.  
 *  
 *  ALPHAI  (output) DOUBLE PRECISION array, dimension (N)  
 *          The imaginary parts of each scalar alpha defining an  
 *          eigenvalue of GNEP.  
 *          If ALPHAI(j) is zero, then the j-th eigenvalue is real; if  
 *          positive, then the j-th and (j+1)-st eigenvalues are a  
 *          complex conjugate pair, with ALPHAI(j+1) = -ALPHAI(j).  
 *  
 *  BETA    (output) DOUBLE PRECISION array, dimension (N)  
 *          The scalars beta that define the eigenvalues of GNEP.  
 *          Together, the quantities alpha = (ALPHAR(j),ALPHAI(j)) and  
 *          beta = BETA(j) represent the j-th eigenvalue of the matrix  
 *          pair (A,B), in one of the forms lambda = alpha/beta or  
 *          mu = beta/alpha.  Since either lambda or mu may overflow,  
 *          they should not, in general, be computed.  
 *  
 *  Q       (input/output) DOUBLE PRECISION array, dimension (LDQ, N)  
 *          On entry, if COMPZ = 'V', the orthogonal matrix Q1 used in  
 *          the reduction of (A,B) to generalized Hessenberg form.  
 *          On exit, if COMPZ = 'I', the orthogonal matrix of left Schur  
 *          vectors of (H,T), and if COMPZ = 'V', the orthogonal matrix  
 *          of left Schur vectors of (A,B).  
 *          Not referenced if COMPZ = 'N'.  
 *  
 *  LDQ     (input) INTEGER  
 *          The leading dimension of the array Q.  LDQ >= 1.  
 *          If COMPQ='V' or 'I', then LDQ >= N.  
 *  
 *  Z       (input/output) DOUBLE PRECISION array, dimension (LDZ, N)  
 *          On entry, if COMPZ = 'V', the orthogonal matrix Z1 used in  
 *          the reduction of (A,B) to generalized Hessenberg form.  
 *          On exit, if COMPZ = 'I', the orthogonal matrix of  
 *          right Schur vectors of (H,T), and if COMPZ = 'V', the  
 *          orthogonal matrix of right Schur vectors of (A,B).  
 *          Not referenced if COMPZ = 'N'.  
 *  
 *  LDZ     (input) INTEGER  
 *          The leading dimension of the array Z.  LDZ >= 1.  
 *          If COMPZ='V' or 'I', then LDZ >= N.  
 *  
 *  WORK    (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LWORK))  
 *          On exit, if INFO >= 0, WORK(1) returns the optimal LWORK.  
 *  
 *  LWORK   (input) INTEGER  
 *          The dimension of the array WORK.  LWORK >= max(1,N).  
 *  
 *          If LWORK = -1, then a workspace query is assumed; the routine  
 *          only calculates the optimal size of the WORK array, returns  
 *          this value as the first entry of the WORK array, and no error  
 *          message related to LWORK is issued by XERBLA.  
 *  
 *  INFO    (output) INTEGER  
 *          = 0: successful exit  
 *          < 0: if INFO = -i, the i-th argument had an illegal value  
 *          = 1,...,N: the QZ iteration did not converge.  (H,T) is not  
 *                     in Schur form, but ALPHAR(i), ALPHAI(i), and  
 *                     BETA(i), i=INFO+1,...,N should be correct.  
 *          = N+1,...,2*N: the shift calculation failed.  (H,T) is not  
 *                     in Schur form, but ALPHAR(i), ALPHAI(i), and  
 *                     BETA(i), i=INFO-N+1,...,N should be correct.  
 *  
 *  Further Details  
 *  ===============  
 *  
 *  Iteration counters:  
 *  
 *  JITER  -- counts iterations.  
 *  IITER  -- counts iterations run since ILAST was last  
 *            changed.  This is therefore reset only when a 1-by-1 or  
 *            2-by-2 block deflates off the bottom.  
 *  
 *  =====================================================================  *  =====================================================================
 *  *
 *     .. Parameters ..  *     .. Parameters ..
Removed from v.1.8  
changed lines
  Added in v.1.9

CVSweb interface <joel.bertrand@systella.fr>