cgels.f.html
来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 448 行 · 第 1/3 页
HTML
448 行
</span><span class="comment">*</span><span class="comment"> .. External Functions ..
</span> LOGICAL <a name="LSAME.133"></a><a href="lsame.f.html#LSAME.1">LSAME</a>
INTEGER <a name="ILAENV.134"></a><a href="hfy-index.html#ILAENV">ILAENV</a>
REAL <a name="CLANGE.135"></a><a href="clange.f.html#CLANGE.1">CLANGE</a>, <a name="SLAMCH.135"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>
EXTERNAL <a name="LSAME.136"></a><a href="lsame.f.html#LSAME.1">LSAME</a>, <a name="ILAENV.136"></a><a href="hfy-index.html#ILAENV">ILAENV</a>, <a name="CLANGE.136"></a><a href="clange.f.html#CLANGE.1">CLANGE</a>, <a name="SLAMCH.136"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. External Subroutines ..
</span> EXTERNAL <a name="CGELQF.139"></a><a href="cgelqf.f.html#CGELQF.1">CGELQF</a>, <a name="CGEQRF.139"></a><a href="cgeqrf.f.html#CGEQRF.1">CGEQRF</a>, <a name="CLASCL.139"></a><a href="clascl.f.html#CLASCL.1">CLASCL</a>, <a name="CLASET.139"></a><a href="claset.f.html#CLASET.1">CLASET</a>, <a name="CTRTRS.139"></a><a href="ctrtrs.f.html#CTRTRS.1">CTRTRS</a>, <a name="CUNMLQ.139"></a><a href="cunmlq.f.html#CUNMLQ.1">CUNMLQ</a>,
$ <a name="CUNMQR.140"></a><a href="cunmqr.f.html#CUNMQR.1">CUNMQR</a>, <a name="SLABAD.140"></a><a href="slabad.f.html#SLABAD.1">SLABAD</a>, <a name="XERBLA.140"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. Intrinsic Functions ..
</span> INTRINSIC MAX, MIN, REAL
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. Executable Statements ..
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Test the input arguments.
</span><span class="comment">*</span><span class="comment">
</span> INFO = 0
MN = MIN( M, N )
LQUERY = ( LWORK.EQ.-1 )
IF( .NOT.( <a name="LSAME.152"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( TRANS, <span class="string">'N'</span> ) .OR. <a name="LSAME.152"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( TRANS, <span class="string">'C'</span> ) ) ) THEN
INFO = -1
ELSE IF( M.LT.0 ) THEN
INFO = -2
ELSE IF( N.LT.0 ) THEN
INFO = -3
ELSE IF( NRHS.LT.0 ) THEN
INFO = -4
ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
INFO = -6
ELSE IF( LDB.LT.MAX( 1, M, N ) ) THEN
INFO = -8
ELSE IF( LWORK.LT.MAX( 1, MN+MAX( MN, NRHS ) ) .AND.
$ .NOT.LQUERY ) THEN
INFO = -10
END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Figure out optimal block size
</span><span class="comment">*</span><span class="comment">
</span> IF( INFO.EQ.0 .OR. INFO.EQ.-10 ) THEN
<span class="comment">*</span><span class="comment">
</span> TPSD = .TRUE.
IF( <a name="LSAME.174"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( TRANS, <span class="string">'N'</span> ) )
$ TPSD = .FALSE.
<span class="comment">*</span><span class="comment">
</span> IF( M.GE.N ) THEN
NB = <a name="ILAENV.178"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="CGEQRF.178"></a><a href="cgeqrf.f.html#CGEQRF.1">CGEQRF</a>'</span>, <span class="string">' '</span>, M, N, -1, -1 )
IF( TPSD ) THEN
NB = MAX( NB, <a name="ILAENV.180"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="CUNMQR.180"></a><a href="cunmqr.f.html#CUNMQR.1">CUNMQR</a>'</span>, <span class="string">'LN'</span>, M, NRHS, N,
$ -1 ) )
ELSE
NB = MAX( NB, <a name="ILAENV.183"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="CUNMQR.183"></a><a href="cunmqr.f.html#CUNMQR.1">CUNMQR</a>'</span>, <span class="string">'LC'</span>, M, NRHS, N,
$ -1 ) )
END IF
ELSE
NB = <a name="ILAENV.187"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="CGELQF.187"></a><a href="cgelqf.f.html#CGELQF.1">CGELQF</a>'</span>, <span class="string">' '</span>, M, N, -1, -1 )
IF( TPSD ) THEN
NB = MAX( NB, <a name="ILAENV.189"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="CUNMLQ.189"></a><a href="cunmlq.f.html#CUNMLQ.1">CUNMLQ</a>'</span>, <span class="string">'LC'</span>, N, NRHS, M,
$ -1 ) )
ELSE
NB = MAX( NB, <a name="ILAENV.192"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="CUNMLQ.192"></a><a href="cunmlq.f.html#CUNMLQ.1">CUNMLQ</a>'</span>, <span class="string">'LN'</span>, N, NRHS, M,
$ -1 ) )
END IF
END IF
<span class="comment">*</span><span class="comment">
</span> WSIZE = MAX( 1, MN + MAX( MN, NRHS )*NB )
WORK( 1 ) = REAL( WSIZE )
<span class="comment">*</span><span class="comment">
</span> END IF
<span class="comment">*</span><span class="comment">
</span> IF( INFO.NE.0 ) THEN
CALL <a name="XERBLA.203"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>( <span class="string">'<a name="CGELS.203"></a><a href="cgels.f.html#CGELS.1">CGELS</a> '</span>, -INFO )
RETURN
ELSE IF( LQUERY ) THEN
RETURN
END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Quick return if possible
</span><span class="comment">*</span><span class="comment">
</span> IF( MIN( M, N, NRHS ).EQ.0 ) THEN
CALL <a name="CLASET.212"></a><a href="claset.f.html#CLASET.1">CLASET</a>( <span class="string">'Full'</span>, MAX( M, N ), NRHS, CZERO, CZERO, B, LDB )
RETURN
END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Get machine parameters
</span><span class="comment">*</span><span class="comment">
</span> SMLNUM = <a name="SLAMCH.218"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>( <span class="string">'S'</span> ) / <a name="SLAMCH.218"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>( <span class="string">'P'</span> )
BIGNUM = ONE / SMLNUM
CALL <a name="SLABAD.220"></a><a href="slabad.f.html#SLABAD.1">SLABAD</a>( SMLNUM, BIGNUM )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Scale A, B if max element outside range [SMLNUM,BIGNUM]
</span><span class="comment">*</span><span class="comment">
</span> ANRM = <a name="CLANGE.224"></a><a href="clange.f.html#CLANGE.1">CLANGE</a>( <span class="string">'M'</span>, M, N, A, LDA, RWORK )
IASCL = 0
IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Scale matrix norm up to SMLNUM
</span><span class="comment">*</span><span class="comment">
</span> CALL <a name="CLASCL.230"></a><a href="clascl.f.html#CLASCL.1">CLASCL</a>( <span class="string">'G'</span>, 0, 0, ANRM, SMLNUM, M, N, A, LDA, INFO )
IASCL = 1
ELSE IF( ANRM.GT.BIGNUM ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Scale matrix norm down to BIGNUM
</span><span class="comment">*</span><span class="comment">
</span> CALL <a name="CLASCL.236"></a><a href="clascl.f.html#CLASCL.1">CLASCL</a>( <span class="string">'G'</span>, 0, 0, ANRM, BIGNUM, M, N, A, LDA, INFO )
IASCL = 2
ELSE IF( ANRM.EQ.ZERO ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Matrix all zero. Return zero solution.
</span><span class="comment">*</span><span class="comment">
</span> CALL <a name="CLASET.242"></a><a href="claset.f.html#CLASET.1">CLASET</a>( <span class="string">'F'</span>, MAX( M, N ), NRHS, CZERO, CZERO, B, LDB )
GO TO 50
END IF
<span class="comment">*</span><span class="comment">
</span> BROW = M
IF( TPSD )
$ BROW = N
BNRM = <a name="CLANGE.249"></a><a href="clange.f.html#CLANGE.1">CLANGE</a>( <span class="string">'M'</span>, BROW, NRHS, B, LDB, RWORK )
IBSCL = 0
IF( BNRM.GT.ZERO .AND. BNRM.LT.SMLNUM ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Scale matrix norm up to SMLNUM
</span><span class="comment">*</span><span class="comment">
</span> CALL <a name="CLASCL.255"></a><a href="clascl.f.html#CLASCL.1">CLASCL</a>( <span class="string">'G'</span>, 0, 0, BNRM, SMLNUM, BROW, NRHS, B, LDB,
$ INFO )
IBSCL = 1
ELSE IF( BNRM.GT.BIGNUM ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Scale matrix norm down to BIGNUM
</span><span class="comment">*</span><span class="comment">
</span> CALL <a name="CLASCL.262"></a><a href="clascl.f.html#CLASCL.1">CLASCL</a>( <span class="string">'G'</span>, 0, 0, BNRM, BIGNUM, BROW, NRHS, B, LDB,
$ INFO )
IBSCL = 2
END IF
<span class="comment">*</span><span class="comment">
</span> IF( M.GE.N ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> compute QR factorization of A
</span><span class="comment">*</span><span class="comment">
</span> CALL <a name="CGEQRF.271"></a><a href="cgeqrf.f.html#CGEQRF.1">CGEQRF</a>( M, N, A, LDA, WORK( 1 ), WORK( MN+1 ), LWORK-MN,
$ INFO )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> workspace at least N, optimally N*NB
</span><span class="comment">*</span><span class="comment">
</span> IF( .NOT.TPSD ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Least-Squares Problem min || A * X - B ||
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> B(1:M,1:NRHS) := Q' * B(1:M,1:NRHS)
</span><span class="comment">*</span><span class="comment">
⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?