Numérotation des lignes :

dqrlss
C DQRLSS    SOURCE    PV090527  23/02/07    11:58:37     11592                subroutine dqrlss ( a, lda, m, n, kr, b, x, r, jpvt, qraux ) c*****************************************************************************80cc! DQRLSS solves a linear system in a least squares sense.cc  Discussion:cc    DQRLSS must be preceeded by a call to DQRANK.cc    The system is to be solved isc      A * X = Bc    wherec      A is an M by N matrix with rank KR, as determined by DQRANK,c      B is a given M-vector,c      X is the N-vector to be computed.cc    A solution X, with at most KR nonzero components, is found whichc    minimizes the 2-norm of the residual (A*X-B).cc    Once the matrix A has been formed, DQRANK should bec    called once to decompose it.  Then, for each right handc    side B, DQRLSS should be called once to obtain thec    solution and residual.cc  Modified:cc    15 April 2012cc  Parameters:cc    Input, real*8 A(LDA,N), the QR factorization informationc    from DQRANK.  The triangular matrix R of the QR factorization isc    contained in the upper triangle and information needed to recoverc    the orthogonal matrix Q is stored below the diagonal in A and inc    the vector QRAUX.cc    Input, integer LDA, the leading dimension of A, which mustc    be at least M.cc    Input, integer M, the number of rows of A.cc    Input, integer N, the number of columns of A.cc    Input, integer KR, the rank of the matrix, as estimatedc    by DQRANK.cc    Input, real*8 B(M), the right hand side of the linear system.cc    Output, real*8 X(N), a least squares solution to thec    linear system.cc    Output, real*8 R(M), the residual, B - A*X.  R mayc    overwite B.cc    Input, integer JPVT(N), the pivot information from DQRANK.c    Columns JPVT(1), ..., JPVT(KR) of the original matrix are linearlyc    independent to within the tolerance TOL and the remaining columnsc    are linearly dependent.cc    Input, real*8 QRAUX(N), auxiliary information from DQRANKc    defining the QR factorization.         implicit real*8 (a-h,o-z)        implicit integer (i-n)         integer lda        integer m        integer n         real*8 a(lda,n)        real*8 b(m)        integer info        integer j        integer job        integer jpvt(n)        integer k        integer kr        real*8 qraux(n)        real*8 r(m)        real*8 t        real*8 x(n)         if ( kr /= 0 ) then          job = 110          call dqrsl ( a, lda, m, kr, qraux, b, r, r, x, r, r, job,info)        end if         jpvt(1:n) = - jpvt(1:n)         x(kr+1:n) = 0.0D+00         do j = 1, n           if ( jpvt(j) .le. 0 ) then             k = -jpvt(j)            jpvt(j) = k             do while ( k /= j )              t = x(j)              x(j) = x(k)              x(k) = t              jpvt(k) = -jpvt(k)              k = jpvt(k)            end do           end if         end do         return      end  

© Cast3M 2003 - Tous droits réservés.
Mentions légales