File:  [local] / rpl / lapack / lapack / dlags2.f
Revision 1.16: download - view: text, annotated - select for diffs - revision graph
Sat Jun 17 10:53:53 2017 UTC (6 years, 10 months ago) by bertrand
Branches: MAIN
CVS tags: HEAD
Mise à jour de lapack.

    1: *> \brief \b DLAGS2 computes 2-by-2 orthogonal matrices U, V, and Q, and applies them to matrices A and B such that the rows of the transformed A and B are parallel.
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at
    6: *            http://www.netlib.org/lapack/explore-html/
    7: *
    8: *> \htmlonly
    9: *> Download DLAGS2 + dependencies
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlags2.f">
   11: *> [TGZ]</a>
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlags2.f">
   13: *> [ZIP]</a>
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlags2.f">
   15: *> [TXT]</a>
   16: *> \endhtmlonly
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE DLAGS2( UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV,
   22: *                          SNV, CSQ, SNQ )
   23: *
   24: *       .. Scalar Arguments ..
   25: *       LOGICAL            UPPER
   26: *       DOUBLE PRECISION   A1, A2, A3, B1, B2, B3, CSQ, CSU, CSV, SNQ,
   27: *      $                   SNU, SNV
   28: *       ..
   29: *
   30: *
   31: *> \par Purpose:
   32: *  =============
   33: *>
   34: *> \verbatim
   35: *>
   36: *> DLAGS2 computes 2-by-2 orthogonal matrices U, V and Q, such
   37: *> that if ( UPPER ) then
   38: *>
   39: *>           U**T *A*Q = U**T *( A1 A2 )*Q = ( x  0  )
   40: *>                             ( 0  A3 )     ( x  x  )
   41: *> and
   42: *>           V**T*B*Q = V**T *( B1 B2 )*Q = ( x  0  )
   43: *>                            ( 0  B3 )     ( x  x  )
   44: *>
   45: *> or if ( .NOT.UPPER ) then
   46: *>
   47: *>           U**T *A*Q = U**T *( A1 0  )*Q = ( x  x  )
   48: *>                             ( A2 A3 )     ( 0  x  )
   49: *> and
   50: *>           V**T*B*Q = V**T*( B1 0  )*Q = ( x  x  )
   51: *>                           ( B2 B3 )     ( 0  x  )
   52: *>
   53: *> The rows of the transformed A and B are parallel, where
   54: *>
   55: *>   U = (  CSU  SNU ), V = (  CSV SNV ), Q = (  CSQ   SNQ )
   56: *>       ( -SNU  CSU )      ( -SNV CSV )      ( -SNQ   CSQ )
   57: *>
   58: *> Z**T denotes the transpose of Z.
   59: *>
   60: *> \endverbatim
   61: *
   62: *  Arguments:
   63: *  ==========
   64: *
   65: *> \param[in] UPPER
   66: *> \verbatim
   67: *>          UPPER is LOGICAL
   68: *>          = .TRUE.: the input matrices A and B are upper triangular.
   69: *>          = .FALSE.: the input matrices A and B are lower triangular.
   70: *> \endverbatim
   71: *>
   72: *> \param[in] A1
   73: *> \verbatim
   74: *>          A1 is DOUBLE PRECISION
   75: *> \endverbatim
   76: *>
   77: *> \param[in] A2
   78: *> \verbatim
   79: *>          A2 is DOUBLE PRECISION
   80: *> \endverbatim
   81: *>
   82: *> \param[in] A3
   83: *> \verbatim
   84: *>          A3 is DOUBLE PRECISION
   85: *>          On entry, A1, A2 and A3 are elements of the input 2-by-2
   86: *>          upper (lower) triangular matrix A.
   87: *> \endverbatim
   88: *>
   89: *> \param[in] B1
   90: *> \verbatim
   91: *>          B1 is DOUBLE PRECISION
   92: *> \endverbatim
   93: *>
   94: *> \param[in] B2
   95: *> \verbatim
   96: *>          B2 is DOUBLE PRECISION
   97: *> \endverbatim
   98: *>
   99: *> \param[in] B3
  100: *> \verbatim
  101: *>          B3 is DOUBLE PRECISION
  102: *>          On entry, B1, B2 and B3 are elements of the input 2-by-2
  103: *>          upper (lower) triangular matrix B.
  104: *> \endverbatim
  105: *>
  106: *> \param[out] CSU
  107: *> \verbatim
  108: *>          CSU is DOUBLE PRECISION
  109: *> \endverbatim
  110: *>
  111: *> \param[out] SNU
  112: *> \verbatim
  113: *>          SNU is DOUBLE PRECISION
  114: *>          The desired orthogonal matrix U.
  115: *> \endverbatim
  116: *>
  117: *> \param[out] CSV
  118: *> \verbatim
  119: *>          CSV is DOUBLE PRECISION
  120: *> \endverbatim
  121: *>
  122: *> \param[out] SNV
  123: *> \verbatim
  124: *>          SNV is DOUBLE PRECISION
  125: *>          The desired orthogonal matrix V.
  126: *> \endverbatim
  127: *>
  128: *> \param[out] CSQ
  129: *> \verbatim
  130: *>          CSQ is DOUBLE PRECISION
  131: *> \endverbatim
  132: *>
  133: *> \param[out] SNQ
  134: *> \verbatim
  135: *>          SNQ is DOUBLE PRECISION
  136: *>          The desired orthogonal matrix Q.
  137: *> \endverbatim
  138: *
  139: *  Authors:
  140: *  ========
  141: *
  142: *> \author Univ. of Tennessee
  143: *> \author Univ. of California Berkeley
  144: *> \author Univ. of Colorado Denver
  145: *> \author NAG Ltd.
  146: *
  147: *> \date December 2016
  148: *
  149: *> \ingroup doubleOTHERauxiliary
  150: *
  151: *  =====================================================================
  152:       SUBROUTINE DLAGS2( UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV,
  153:      $                   SNV, CSQ, SNQ )
  154: *
  155: *  -- LAPACK auxiliary routine (version 3.7.0) --
  156: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  157: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  158: *     December 2016
  159: *
  160: *     .. Scalar Arguments ..
  161:       LOGICAL            UPPER
  162:       DOUBLE PRECISION   A1, A2, A3, B1, B2, B3, CSQ, CSU, CSV, SNQ,
  163:      $                   SNU, SNV
  164: *     ..
  165: *
  166: *  =====================================================================
  167: *
  168: *     .. Parameters ..
  169:       DOUBLE PRECISION   ZERO
  170:       PARAMETER          ( ZERO = 0.0D+0 )
  171: *     ..
  172: *     .. Local Scalars ..
  173:       DOUBLE PRECISION   A, AUA11, AUA12, AUA21, AUA22, AVB11, AVB12,
  174:      $                   AVB21, AVB22, B, C, CSL, CSR, D, R, S1, S2,
  175:      $                   SNL, SNR, UA11, UA11R, UA12, UA21, UA22, UA22R,
  176:      $                   VB11, VB11R, VB12, VB21, VB22, VB22R
  177: *     ..
  178: *     .. External Subroutines ..
  179:       EXTERNAL           DLARTG, DLASV2
  180: *     ..
  181: *     .. Intrinsic Functions ..
  182:       INTRINSIC          ABS
  183: *     ..
  184: *     .. Executable Statements ..
  185: *
  186:       IF( UPPER ) THEN
  187: *
  188: *        Input matrices A and B are upper triangular matrices
  189: *
  190: *        Form matrix C = A*adj(B) = ( a b )
  191: *                                   ( 0 d )
  192: *
  193:          A = A1*B3
  194:          D = A3*B1
  195:          B = A2*B1 - A1*B2
  196: *
  197: *        The SVD of real 2-by-2 triangular C
  198: *
  199: *         ( CSL -SNL )*( A B )*(  CSR  SNR ) = ( R 0 )
  200: *         ( SNL  CSL ) ( 0 D ) ( -SNR  CSR )   ( 0 T )
  201: *
  202:          CALL DLASV2( A, B, D, S1, S2, SNR, CSR, SNL, CSL )
  203: *
  204:          IF( ABS( CSL ).GE.ABS( SNL ) .OR. ABS( CSR ).GE.ABS( SNR ) )
  205:      $        THEN
  206: *
  207: *           Compute the (1,1) and (1,2) elements of U**T *A and V**T *B,
  208: *           and (1,2) element of |U|**T *|A| and |V|**T *|B|.
  209: *
  210:             UA11R = CSL*A1
  211:             UA12 = CSL*A2 + SNL*A3
  212: *
  213:             VB11R = CSR*B1
  214:             VB12 = CSR*B2 + SNR*B3
  215: *
  216:             AUA12 = ABS( CSL )*ABS( A2 ) + ABS( SNL )*ABS( A3 )
  217:             AVB12 = ABS( CSR )*ABS( B2 ) + ABS( SNR )*ABS( B3 )
  218: *
  219: *           zero (1,2) elements of U**T *A and V**T *B
  220: *
  221:             IF( ( ABS( UA11R )+ABS( UA12 ) ).NE.ZERO ) THEN
  222:                IF( AUA12 / ( ABS( UA11R )+ABS( UA12 ) ).LE.AVB12 /
  223:      $             ( ABS( VB11R )+ABS( VB12 ) ) ) THEN
  224:                   CALL DLARTG( -UA11R, UA12, CSQ, SNQ, R )
  225:                ELSE
  226:                   CALL DLARTG( -VB11R, VB12, CSQ, SNQ, R )
  227:                END IF
  228:             ELSE
  229:                CALL DLARTG( -VB11R, VB12, CSQ, SNQ, R )
  230:             END IF
  231: *
  232:             CSU = CSL
  233:             SNU = -SNL
  234:             CSV = CSR
  235:             SNV = -SNR
  236: *
  237:          ELSE
  238: *
  239: *           Compute the (2,1) and (2,2) elements of U**T *A and V**T *B,
  240: *           and (2,2) element of |U|**T *|A| and |V|**T *|B|.
  241: *
  242:             UA21 = -SNL*A1
  243:             UA22 = -SNL*A2 + CSL*A3
  244: *
  245:             VB21 = -SNR*B1
  246:             VB22 = -SNR*B2 + CSR*B3
  247: *
  248:             AUA22 = ABS( SNL )*ABS( A2 ) + ABS( CSL )*ABS( A3 )
  249:             AVB22 = ABS( SNR )*ABS( B2 ) + ABS( CSR )*ABS( B3 )
  250: *
  251: *           zero (2,2) elements of U**T*A and V**T*B, and then swap.
  252: *
  253:             IF( ( ABS( UA21 )+ABS( UA22 ) ).NE.ZERO ) THEN
  254:                IF( AUA22 / ( ABS( UA21 )+ABS( UA22 ) ).LE.AVB22 /
  255:      $             ( ABS( VB21 )+ABS( VB22 ) ) ) THEN
  256:                   CALL DLARTG( -UA21, UA22, CSQ, SNQ, R )
  257:                ELSE
  258:                   CALL DLARTG( -VB21, VB22, CSQ, SNQ, R )
  259:                END IF
  260:             ELSE
  261:                CALL DLARTG( -VB21, VB22, CSQ, SNQ, R )
  262:             END IF
  263: *
  264:             CSU = SNL
  265:             SNU = CSL
  266:             CSV = SNR
  267:             SNV = CSR
  268: *
  269:          END IF
  270: *
  271:       ELSE
  272: *
  273: *        Input matrices A and B are lower triangular matrices
  274: *
  275: *        Form matrix C = A*adj(B) = ( a 0 )
  276: *                                   ( c d )
  277: *
  278:          A = A1*B3
  279:          D = A3*B1
  280:          C = A2*B3 - A3*B2
  281: *
  282: *        The SVD of real 2-by-2 triangular C
  283: *
  284: *         ( CSL -SNL )*( A 0 )*(  CSR  SNR ) = ( R 0 )
  285: *         ( SNL  CSL ) ( C D ) ( -SNR  CSR )   ( 0 T )
  286: *
  287:          CALL DLASV2( A, C, D, S1, S2, SNR, CSR, SNL, CSL )
  288: *
  289:          IF( ABS( CSR ).GE.ABS( SNR ) .OR. ABS( CSL ).GE.ABS( SNL ) )
  290:      $        THEN
  291: *
  292: *           Compute the (2,1) and (2,2) elements of U**T *A and V**T *B,
  293: *           and (2,1) element of |U|**T *|A| and |V|**T *|B|.
  294: *
  295:             UA21 = -SNR*A1 + CSR*A2
  296:             UA22R = CSR*A3
  297: *
  298:             VB21 = -SNL*B1 + CSL*B2
  299:             VB22R = CSL*B3
  300: *
  301:             AUA21 = ABS( SNR )*ABS( A1 ) + ABS( CSR )*ABS( A2 )
  302:             AVB21 = ABS( SNL )*ABS( B1 ) + ABS( CSL )*ABS( B2 )
  303: *
  304: *           zero (2,1) elements of U**T *A and V**T *B.
  305: *
  306:             IF( ( ABS( UA21 )+ABS( UA22R ) ).NE.ZERO ) THEN
  307:                IF( AUA21 / ( ABS( UA21 )+ABS( UA22R ) ).LE.AVB21 /
  308:      $             ( ABS( VB21 )+ABS( VB22R ) ) ) THEN
  309:                   CALL DLARTG( UA22R, UA21, CSQ, SNQ, R )
  310:                ELSE
  311:                   CALL DLARTG( VB22R, VB21, CSQ, SNQ, R )
  312:                END IF
  313:             ELSE
  314:                CALL DLARTG( VB22R, VB21, CSQ, SNQ, R )
  315:             END IF
  316: *
  317:             CSU = CSR
  318:             SNU = -SNR
  319:             CSV = CSL
  320:             SNV = -SNL
  321: *
  322:          ELSE
  323: *
  324: *           Compute the (1,1) and (1,2) elements of U**T *A and V**T *B,
  325: *           and (1,1) element of |U|**T *|A| and |V|**T *|B|.
  326: *
  327:             UA11 = CSR*A1 + SNR*A2
  328:             UA12 = SNR*A3
  329: *
  330:             VB11 = CSL*B1 + SNL*B2
  331:             VB12 = SNL*B3
  332: *
  333:             AUA11 = ABS( CSR )*ABS( A1 ) + ABS( SNR )*ABS( A2 )
  334:             AVB11 = ABS( CSL )*ABS( B1 ) + ABS( SNL )*ABS( B2 )
  335: *
  336: *           zero (1,1) elements of U**T*A and V**T*B, and then swap.
  337: *
  338:             IF( ( ABS( UA11 )+ABS( UA12 ) ).NE.ZERO ) THEN
  339:                IF( AUA11 / ( ABS( UA11 )+ABS( UA12 ) ).LE.AVB11 /
  340:      $             ( ABS( VB11 )+ABS( VB12 ) ) ) THEN
  341:                   CALL DLARTG( UA12, UA11, CSQ, SNQ, R )
  342:                ELSE
  343:                   CALL DLARTG( VB12, VB11, CSQ, SNQ, R )
  344:                END IF
  345:             ELSE
  346:                CALL DLARTG( VB12, VB11, CSQ, SNQ, R )
  347:             END IF
  348: *
  349:             CSU = SNR
  350:             SNU = CSR
  351:             CSV = SNL
  352:             SNV = CSL
  353: *
  354:          END IF
  355: *
  356:       END IF
  357: *
  358:       RETURN
  359: *
  360: *     End of DLAGS2
  361: *
  362:       END

CVSweb interface <joel.bertrand@systella.fr>