ctgsyl.f.html

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

HTML
597 行
字号
123
      GO TO 60
<span class="comment">*</span><span class="comment">
</span>   70 CONTINUE
      IWORK( Q+1 ) = N + 1
      IF( IWORK( Q ).EQ.IWORK( Q+1 ) )
     $   Q = Q - 1
<span class="comment">*</span><span class="comment">
</span>      IF( NOTRAN ) THEN
         DO 150 IROUND = 1, ISOLVE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Solve (I, J) - subsystem
</span><span class="comment">*</span><span class="comment">               A(I, I) * R(I, J) - L(I, J) * B(J, J) = C(I, J)
</span><span class="comment">*</span><span class="comment">               D(I, I) * R(I, J) - L(I, J) * E(J, J) = F(I, J)
</span><span class="comment">*</span><span class="comment">           for I = P, P - 1, ..., 1; J = 1, 2, ..., Q
</span><span class="comment">*</span><span class="comment">
</span>            PQ = 0
            SCALE = ONE
            DSCALE = ZERO
            DSUM = ONE
            DO 130 J = P + 2, Q
               JS = IWORK( J )
               JE = IWORK( J+1 ) - 1
               NB = JE - JS + 1
               DO 120 I = P, 1, -1
                  IS = IWORK( I )
                  IE = IWORK( I+1 ) - 1
                  MB = IE - IS + 1
                  CALL <a name="CTGSY2.407"></a><a href="ctgsy2.f.html#CTGSY2.1">CTGSY2</a>( TRANS, IFUNC, MB, NB, A( IS, IS ), LDA,
     $                         B( JS, JS ), LDB, C( IS, JS ), LDC,
     $                         D( IS, IS ), LDD, E( JS, JS ), LDE,
     $                         F( IS, JS ), LDF, SCALOC, DSUM, DSCALE,
     $                         LINFO )
                  IF( LINFO.GT.0 )
     $               INFO = LINFO
                  PQ = PQ + MB*NB
                  IF( SCALOC.NE.ONE ) THEN
                     DO 80 K = 1, JS - 1
                        CALL CSCAL( M, CMPLX( SCALOC, ZERO ), C( 1, K ),
     $                              1 )
                        CALL CSCAL( M, CMPLX( SCALOC, ZERO ), F( 1, K ),
     $                              1 )
   80                CONTINUE
                     DO 90 K = JS, JE
                        CALL CSCAL( IS-1, CMPLX( SCALOC, ZERO ),
     $                              C( 1, K ), 1 )
                        CALL CSCAL( IS-1, CMPLX( SCALOC, ZERO ),
     $                              F( 1, K ), 1 )
   90                CONTINUE
                     DO 100 K = JS, JE
                        CALL CSCAL( M-IE, CMPLX( SCALOC, ZERO ),
     $                              C( IE+1, K ), 1 )
                        CALL CSCAL( M-IE, CMPLX( SCALOC, ZERO ),
     $                              F( IE+1, K ), 1 )
  100                CONTINUE
                     DO 110 K = JE + 1, N
                        CALL CSCAL( M, CMPLX( SCALOC, ZERO ), C( 1, K ),
     $                              1 )
                        CALL CSCAL( M, CMPLX( SCALOC, ZERO ), F( 1, K ),
     $                              1 )
  110                CONTINUE
                     SCALE = SCALE*SCALOC
                  END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 Substitute R(I,J) and L(I,J) into remaining equation.
</span><span class="comment">*</span><span class="comment">
</span>                  IF( I.GT.1 ) THEN
                     CALL CGEMM( <span class="string">'N'</span>, <span class="string">'N'</span>, IS-1, NB, MB,
     $                           CMPLX( -ONE, ZERO ), A( 1, IS ), LDA,
     $                           C( IS, JS ), LDC, CMPLX( ONE, ZERO ),
     $                           C( 1, JS ), LDC )
                     CALL CGEMM( <span class="string">'N'</span>, <span class="string">'N'</span>, IS-1, NB, MB,
     $                           CMPLX( -ONE, ZERO ), D( 1, IS ), LDD,
     $                           C( IS, JS ), LDC, CMPLX( ONE, ZERO ),
     $                           F( 1, JS ), LDF )
                  END IF
                  IF( J.LT.Q ) THEN
                     CALL CGEMM( <span class="string">'N'</span>, <span class="string">'N'</span>, MB, N-JE, NB,
     $                           CMPLX( ONE, ZERO ), F( IS, JS ), LDF,
     $                           B( JS, JE+1 ), LDB, CMPLX( ONE, ZERO ),
     $                           C( IS, JE+1 ), LDC )
                     CALL CGEMM( <span class="string">'N'</span>, <span class="string">'N'</span>, MB, N-JE, NB,
     $                           CMPLX( ONE, ZERO ), F( IS, JS ), LDF,
     $                           E( JS, JE+1 ), LDE, CMPLX( ONE, ZERO ),
     $                           F( IS, JE+1 ), LDF )
                  END IF
  120          CONTINUE
  130       CONTINUE
            IF( DSCALE.NE.ZERO ) THEN
               IF( IJOB.EQ.1 .OR. IJOB.EQ.3 ) THEN
                  DIF = SQRT( REAL( 2*M*N ) ) / ( DSCALE*SQRT( DSUM ) )
               ELSE
                  DIF = SQRT( REAL( PQ ) ) / ( DSCALE*SQRT( DSUM ) )
               END IF
            END IF
            IF( ISOLVE.EQ.2 .AND. IROUND.EQ.1 ) THEN
               IF( NOTRAN ) THEN
                  IFUNC = IJOB
               END IF
               SCALE2 = SCALE
               CALL <a name="CLACPY.479"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'F'</span>, M, N, C, LDC, WORK, M )
               CALL <a name="CLACPY.480"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'F'</span>, M, N, F, LDF, WORK( M*N+1 ), M )
               CALL <a name="CLASET.481"></a><a href="claset.f.html#CLASET.1">CLASET</a>( <span class="string">'F'</span>, M, N, CZERO, CZERO, C, LDC )
               CALL <a name="CLASET.482"></a><a href="claset.f.html#CLASET.1">CLASET</a>( <span class="string">'F'</span>, M, N, CZERO, CZERO, F, LDF )
            ELSE IF( ISOLVE.EQ.2 .AND. IROUND.EQ.2 ) THEN
               CALL <a name="CLACPY.484"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'F'</span>, M, N, WORK, M, C, LDC )
               CALL <a name="CLACPY.485"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'F'</span>, M, N, WORK( M*N+1 ), M, F, LDF )
               SCALE = SCALE2
            END IF
  150    CONTINUE
      ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Solve transposed (I, J)-subsystem
</span><span class="comment">*</span><span class="comment">            A(I, I)' * R(I, J) + D(I, I)' * L(I, J) = C(I, J)
</span><span class="comment">*</span><span class="comment">            R(I, J) * B(J, J)  + L(I, J) * E(J, J) = -F(I, J)
</span><span class="comment">*</span><span class="comment">        for I = 1,2,..., P; J = Q, Q-1,..., 1
</span><span class="comment">*</span><span class="comment">
</span>         SCALE = ONE
         DO 210 I = 1, P
            IS = IWORK( I )
            IE = IWORK( I+1 ) - 1
            MB = IE - IS + 1
            DO 200 J = Q, P + 2, -1
               JS = IWORK( J )
               JE = IWORK( J+1 ) - 1
               NB = JE - JS + 1
               CALL <a name="CTGSY2.505"></a><a href="ctgsy2.f.html#CTGSY2.1">CTGSY2</a>( TRANS, IFUNC, MB, NB, A( IS, IS ), LDA,
     $                      B( JS, JS ), LDB, C( IS, JS ), LDC,
     $                      D( IS, IS ), LDD, E( JS, JS ), LDE,
     $                      F( IS, JS ), LDF, SCALOC, DSUM, DSCALE,
     $                      LINFO )
               IF( LINFO.GT.0 )
     $            INFO = LINFO
               IF( SCALOC.NE.ONE ) THEN
                  DO 160 K = 1, JS - 1
                     CALL CSCAL( M, CMPLX( SCALOC, ZERO ), C( 1, K ),
     $                           1 )
                     CALL CSCAL( M, CMPLX( SCALOC, ZERO ), F( 1, K ),
     $                           1 )
  160             CONTINUE
                  DO 170 K = JS, JE
                     CALL CSCAL( IS-1, CMPLX( SCALOC, ZERO ), C( 1, K ),
     $                           1 )
                     CALL CSCAL( IS-1, CMPLX( SCALOC, ZERO ), F( 1, K ),
     $                           1 )
  170             CONTINUE
                  DO 180 K = JS, JE
                     CALL CSCAL( M-IE, CMPLX( SCALOC, ZERO ),
     $                           C( IE+1, K ), 1 )
                     CALL CSCAL( M-IE, CMPLX( SCALOC, ZERO ),
     $                           F( IE+1, K ), 1 )
  180             CONTINUE
                  DO 190 K = JE + 1, N
                     CALL CSCAL( M, CMPLX( SCALOC, ZERO ), C( 1, K ),
     $                           1 )
                     CALL CSCAL( M, CMPLX( SCALOC, ZERO ), F( 1, K ),
     $                           1 )
  190             CONTINUE
                  SCALE = SCALE*SCALOC
               END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Substitute R(I,J) and L(I,J) into remaining equation.
</span><span class="comment">*</span><span class="comment">
</span>               IF( J.GT.P+2 ) THEN
                  CALL CGEMM( <span class="string">'N'</span>, <span class="string">'C'</span>, MB, JS-1, NB,
     $                        CMPLX( ONE, ZERO ), C( IS, JS ), LDC,
     $                        B( 1, JS ), LDB, CMPLX( ONE, ZERO ),
     $                        F( IS, 1 ), LDF )
                  CALL CGEMM( <span class="string">'N'</span>, <span class="string">'C'</span>, MB, JS-1, NB,
     $                        CMPLX( ONE, ZERO ), F( IS, JS ), LDF,
     $                        E( 1, JS ), LDE, CMPLX( ONE, ZERO ),
     $                        F( IS, 1 ), LDF )
               END IF
               IF( I.LT.P ) THEN
                  CALL CGEMM( <span class="string">'C'</span>, <span class="string">'N'</span>, M-IE, NB, MB,
     $                        CMPLX( -ONE, ZERO ), A( IS, IE+1 ), LDA,
     $                        C( IS, JS ), LDC, CMPLX( ONE, ZERO ),
     $                        C( IE+1, JS ), LDC )
                  CALL CGEMM( <span class="string">'C'</span>, <span class="string">'N'</span>, M-IE, NB, MB,
     $                        CMPLX( -ONE, ZERO ), D( IS, IE+1 ), LDD,
     $                        F( IS, JS ), LDF, CMPLX( ONE, ZERO ),
     $                        C( IE+1, JS ), LDC )
               END IF
  200       CONTINUE
  210    CONTINUE
      END IF
<span class="comment">*</span><span class="comment">
</span>      WORK( 1 ) = LWMIN
<span class="comment">*</span><span class="comment">
</span>      RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     End of <a name="CTGSYL.570"></a><a href="ctgsyl.f.html#CTGSYL.1">CTGSYL</a>
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>
123

ctgsyl.f.html - 源码说明

本页面展示了「famous linear algebra library (LAPACK) ports to windows」中的 ctgsyl.f.html 源码文件，采用 HTML 编程语言编写，共 597 行代码。您可以在线阅读完整代码内容，也可以返回资源详情页下载完整源码包进行本地学习和开发。

虫虫下载站收录了大量与algebra相关的技术资源，包括源代码、技术文档、电路图等，是电子工程师和嵌入式开发者的专业学习平台。

⌨️ 快捷键说明

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