zpbtrf.f.html

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

                  INFO = I + II - 1
                  GO TO 150
               END IF
               IF( I+IB.LE.N ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 Update the relevant part of the trailing submatrix.
</span><span class="comment">*</span><span class="comment">                 If A11 denotes the diagonal block which has just been
</span><span class="comment">*</span><span class="comment">                 factorized, then we need to update the remaining
</span><span class="comment">*</span><span class="comment">                 blocks in the diagram:
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    A11   A12   A13
</span><span class="comment">*</span><span class="comment">                          A22   A23
</span><span class="comment">*</span><span class="comment">                                A33
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 The numbers of rows and columns in the partitioning
</span><span class="comment">*</span><span class="comment">                 are IB, I2, I3 respectively. The blocks A12, A22 and
</span><span class="comment">*</span><span class="comment">                 A23 are empty if IB = KD. The upper triangle of A13
</span><span class="comment">*</span><span class="comment">                 lies outside the band.
</span><span class="comment">*</span><span class="comment">
</span>                  I2 = MIN( KD-IB, N-I-IB+1 )
                  I3 = MIN( IB, N-I-KD+1 )
<span class="comment">*</span><span class="comment">
</span>                  IF( I2.GT.0 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A12
</span><span class="comment">*</span><span class="comment">
</span>                     CALL ZTRSM( <span class="string">'Left'</span>, <span class="string">'Upper'</span>, <span class="string">'Conjugate transpose'</span>,
     $                           <span class="string">'Non-unit'</span>, IB, I2, CONE,
     $                           AB( KD+1, I ), LDAB-1,
     $                           AB( KD+1-IB, I+IB ), LDAB-1 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A22
</span><span class="comment">*</span><span class="comment">
</span>                     CALL ZHERK( <span class="string">'Upper'</span>, <span class="string">'Conjugate transpose'</span>, I2, IB,
     $                           -ONE, AB( KD+1-IB, I+IB ), LDAB-1, ONE,
     $                           AB( KD+1, I+IB ), LDAB-1 )
                  END IF
<span class="comment">*</span><span class="comment">
</span>                  IF( I3.GT.0 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Copy the lower triangle of A13 into the work array.
</span><span class="comment">*</span><span class="comment">
</span>                     DO 40 JJ = 1, I3
                        DO 30 II = JJ, IB
                           WORK( II, JJ ) = AB( II-JJ+1, JJ+I+KD-1 )
   30                   CONTINUE
   40                CONTINUE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A13 (in the work array).
</span><span class="comment">*</span><span class="comment">
</span>                     CALL ZTRSM( <span class="string">'Left'</span>, <span class="string">'Upper'</span>, <span class="string">'Conjugate transpose'</span>,
     $                           <span class="string">'Non-unit'</span>, IB, I3, CONE,
     $                           AB( KD+1, I ), LDAB-1, WORK, LDWORK )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A23
</span><span class="comment">*</span><span class="comment">
</span>                     IF( I2.GT.0 )
     $                  CALL ZGEMM( <span class="string">'Conjugate transpose'</span>,
     $                              <span class="string">'No transpose'</span>, I2, I3, IB, -CONE,
     $                              AB( KD+1-IB, I+IB ), LDAB-1, WORK,
     $                              LDWORK, CONE, AB( 1+IB, I+KD ),
     $                              LDAB-1 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A33
</span><span class="comment">*</span><span class="comment">
</span>                     CALL ZHERK( <span class="string">'Upper'</span>, <span class="string">'Conjugate transpose'</span>, I3, IB,
     $                           -ONE, WORK, LDWORK, ONE,
     $                           AB( KD+1, I+KD ), LDAB-1 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Copy the lower triangle of A13 back into place.
</span><span class="comment">*</span><span class="comment">
</span>                     DO 60 JJ = 1, I3
                        DO 50 II = JJ, IB
                           AB( II-JJ+1, JJ+I+KD-1 ) = WORK( II, JJ )
   50                   CONTINUE
   60                CONTINUE
                  END IF
               END IF
   70       CONTINUE
         ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Compute the Cholesky factorization of a Hermitian band
</span><span class="comment">*</span><span class="comment">           matrix, given the lower triangle of the matrix in band
</span><span class="comment">*</span><span class="comment">           storage.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Zero the lower triangle of the work array.
</span><span class="comment">*</span><span class="comment">
</span>            DO 90 J = 1, NB
               DO 80 I = J + 1, NB
                  WORK( I, J ) = ZERO
   80          CONTINUE
   90       CONTINUE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Process the band matrix one diagonal block at a time.
</span><span class="comment">*</span><span class="comment">
</span>            DO 140 I = 1, N, NB
               IB = MIN( NB, N-I+1 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Factorize the diagonal block
</span><span class="comment">*</span><span class="comment">
</span>               CALL <a name="ZPOTF2.280"></a><a href="zpotf2.f.html#ZPOTF2.1">ZPOTF2</a>( UPLO, IB, AB( 1, I ), LDAB-1, II )
               IF( II.NE.0 ) THEN
                  INFO = I + II - 1
                  GO TO 150
               END IF
               IF( I+IB.LE.N ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 Update the relevant part of the trailing submatrix.
</span><span class="comment">*</span><span class="comment">                 If A11 denotes the diagonal block which has just been
</span><span class="comment">*</span><span class="comment">                 factorized, then we need to update the remaining
</span><span class="comment">*</span><span class="comment">                 blocks in the diagram:
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    A11
</span><span class="comment">*</span><span class="comment">                    A21   A22
</span><span class="comment">*</span><span class="comment">                    A31   A32   A33
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 The numbers of rows and columns in the partitioning
</span><span class="comment">*</span><span class="comment">                 are IB, I2, I3 respectively. The blocks A21, A22 and
</span><span class="comment">*</span><span class="comment">                 A32 are empty if IB = KD. The lower triangle of A31
</span><span class="comment">*</span><span class="comment">                 lies outside the band.
</span><span class="comment">*</span><span class="comment">
</span>                  I2 = MIN( KD-IB, N-I-IB+1 )
                  I3 = MIN( IB, N-I-KD+1 )
<span class="comment">*</span><span class="comment">
</span>                  IF( I2.GT.0 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A21
</span><span class="comment">*</span><span class="comment">
</span>                     CALL ZTRSM( <span class="string">'Right'</span>, <span class="string">'Lower'</span>,
     $                           <span class="string">'Conjugate transpose'</span>, <span class="string">'Non-unit'</span>, I2,
     $                           IB, CONE, AB( 1, I ), LDAB-1,
     $                           AB( 1+IB, I ), LDAB-1 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A22
</span><span class="comment">*</span><span class="comment">
</span>                     CALL ZHERK( <span class="string">'Lower'</span>, <span class="string">'No transpose'</span>, I2, IB, -ONE,
     $                           AB( 1+IB, I ), LDAB-1, ONE,
     $                           AB( 1, I+IB ), LDAB-1 )
                  END IF
<span class="comment">*</span><span class="comment">
</span>                  IF( I3.GT.0 ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Copy the upper triangle of A31 into the work array.
</span><span class="comment">*</span><span class="comment">
</span>                     DO 110 JJ = 1, IB
                        DO 100 II = 1, MIN( JJ, I3 )
                           WORK( II, JJ ) = AB( KD+1-JJ+II, JJ+I-1 )
  100                   CONTINUE
  110                CONTINUE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A31 (in the work array).
</span><span class="comment">*</span><span class="comment">
</span>                     CALL ZTRSM( <span class="string">'Right'</span>, <span class="string">'Lower'</span>,
     $                           <span class="string">'Conjugate transpose'</span>, <span class="string">'Non-unit'</span>, I3,
     $                           IB, CONE, AB( 1, I ), LDAB-1, WORK,
     $                           LDWORK )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A32
</span><span class="comment">*</span><span class="comment">
</span>                     IF( I2.GT.0 )
     $                  CALL ZGEMM( <span class="string">'No transpose'</span>,
     $                              <span class="string">'Conjugate transpose'</span>, I3, I2, IB,
     $                              -CONE, WORK, LDWORK, AB( 1+IB, I ),
     $                              LDAB-1, CONE, AB( 1+KD-IB, I+IB ),
     $                              LDAB-1 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Update A33
</span><span class="comment">*</span><span class="comment">
</span>                     CALL ZHERK( <span class="string">'Lower'</span>, <span class="string">'No transpose'</span>, I3, IB, -ONE,
     $                           WORK, LDWORK, ONE, AB( 1, I+KD ),
     $                           LDAB-1 )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                    Copy the upper triangle of A31 back into place.
</span><span class="comment">*</span><span class="comment">
</span>                     DO 130 JJ = 1, IB
                        DO 120 II = 1, MIN( JJ, I3 )
                           AB( KD+1-JJ+II, JJ+I-1 ) = WORK( II, JJ )
  120                   CONTINUE
  130                CONTINUE
                  END IF
               END IF
  140       CONTINUE
         END IF
      END IF
      RETURN
<span class="comment">*</span><span class="comment">
</span>  150 CONTINUE
      RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     End of <a name="ZPBTRF.369"></a><a href="zpbtrf.f.html#ZPBTRF.1">ZPBTRF</a>
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>