lu7rnk Subroutine

private subroutine lu7rnk(m, n, jsing, lena, parmlu, lenL, lenU, lrow, nrank, a, indc, indr, p, q, lenr, locc, locr, inform, diag)

Arguments

Type IntentOptional Attributes Name
integer(kind=ip), intent(in) :: m
integer(kind=ip), intent(in) :: n
integer(kind=ip), intent(in) :: jsing
integer(kind=ip), intent(in) :: lena
real(kind=rp), intent(inout) :: parmlu(30)
integer(kind=ip), intent(inout) :: lenL
integer(kind=ip), intent(inout) :: lenU
integer(kind=ip), intent(inout) :: lrow
integer(kind=ip), intent(inout) :: nrank
real(kind=rp), intent(inout) :: a(lena)
integer(kind=ip), intent(inout) :: indc(lena)
integer(kind=ip), intent(inout) :: indr(lena)
integer(kind=ip), intent(in) :: p(m)
integer(kind=ip), intent(inout) :: q(n)
integer(kind=ip), intent(inout) :: lenr(m)
integer(kind=ip), intent(inout) :: locc(n)
integer(kind=ip), intent(inout) :: locr(m)
integer(kind=ip), intent(out) :: inform
real(kind=rp), intent(out) :: diag

Called by

proc~~lu7rnk~~CalledByGraph proc~lu7rnk lusol::lu7rnk proc~lu8rpc lusol::lu8rpc proc~lu8rpc->proc~lu7rnk

Source Code

  subroutine lu7rnk( m, n, jsing, lena, parmlu,             &
                     lenL, lenU, lrow, nrank,               &
                     a, indc, indr, p, q, lenr, locc, locr, &
                     inform, diag )

    integer(ip),   intent(in)    :: m, n, jsing, lena,      &
                                    p(m)
    integer(ip),   intent(inout) :: lenL, lenU, lrow, nrank, &
                                    indc(lena), indr(lena), q(n),        &
                                    lenr(m), locc(n), locr(m)
    real(rp),      intent(inout) :: parmlu(30)  ! not used
    real(rp),      intent(inout) :: a(lena)
    integer(ip),   intent(out)   :: inform
    real(rp),      intent(out)   :: diag

    !------------------------------------------------------------------
    ! lu7rnk (check rank) assumes U is currently nrank by n
    ! and determines if row nrank contains an acceptable pivot.
    ! If not, the row is deleted and nrank is decreased by 1.

    ! jsing is an input parameter (not altered).  If jsing is positive,
    ! column jsing has already been judged dependent.  A substitute
    ! (if any) must be some other column.
    !
    ! On exit,
    ! inform = -1 if nrank decreases by 1
    !        =  0 if nrank stays the same
    !        =  1 if there's a fatal error.  (Currently we stop.)
    !
    ! -- Jul 1987: First version.
    ! 09 May 1988: First f77 version.
    ! 13 Dec 2011: First f90 version.
    ! 01 Jan 2012: luparm not used.
    !------------------------------------------------------------------

    integer(ip)             :: iw, jmax, kmax, l, l1, l2, lenw, lmax
    real(rp)                :: Umax, Utol1

    Utol1    = parmlu(4)
    diag     = zero

    ! Find Umax, the largest element in row nrank.

    iw       = p(nrank)
    lenw     = lenr(iw)
    if (lenw == 0) go to 400
    l1       = locr(iw)
    l2       = l1 + lenw - 1
    Umax     = zero
    lmax     = l1

    do l = l1, l2
       if (Umax < abs(a(l))) then
          Umax  = abs(a(l))
          lmax  = l
       end if
    end do

    ! Find which column that guy is in (in pivotal order).
    ! Interchange him with column nrank, then move him to be
    ! the new diagonal at the front of row nrank.

    diag   = a(lmax)
    jmax   = indr(lmax)
    l      = 0

    do kmax = nrank, n
       if (q(kmax) == jmax) then
          l = kmax   ! l allows check below for fatal error
          exit
       end if
    end do

    if (l == 0) go to 800   ! Fatal error

    q(kmax)    = q(nrank)
    q(nrank)   = jmax
    a(lmax)    = a(l1)
    a(l1)      = diag
    indr(lmax) = indr(l1)
    indr(l1)   = jmax

    ! See if the new diagonal is big enough.

    if (Umax <= Utol1) go to 400
    if (jmax == jsing) go to 400

    !------------------------------------------------------------------
    ! The rank stays the same.
    !------------------------------------------------------------------
    inform = 0
    go to 900

    !------------------------------------------------------------------
    ! The rank decreases by one.
    !------------------------------------------------------------------
400 inform = -1
    nrank  = nrank - 1

    if (lenw > 0) then       ! Delete row nrank from U.
       lenU     = lenU - lenw
       lenr(iw) = 0
       do l = l1, l2
          indr(l) = 0
       end do

       if (l2 == lrow) then
          ! This row was at the end of the data structure.
          ! We have to reset lrow.
          ! Preceding rows might already have been deleted, so we
          ! have to be prepared to go all the way back to 1.

          do l = 1, l2
             if (indr(lrow) > 0) go to 900
             lrow  = lrow - 1
          end do
       end if
    end if
    go to 900

    ! 15 Dec 2011: Fatal error (should never happen!).
    ! This is a safeguard during work on the f90 version.

800 inform = 1
    write(*,*) 'Fatal error in LUSOL lu7rnk.  Stopping now'
    stop

900 return

  end subroutine lu7rnk