balance Subroutine

  1. eigensystem.f90
  2. eigensystem
  3. balance

public subroutine balance(n, mat, scal, low, high)

Balances the matrix so that the rows with zero entries off the diagonal are isolated and the remaining columns and rows are resized to have one norm close to 1.

Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: n

Dimension of mat
real(kind=double) :: mat(0:n-1,0:n-1)

nxn scaled matrix
real(kind=double) :: scal(0:n-1)

Contains isolated eigenvalue in the positions 0-low and high-n-1 its other components contain the scaling factors for transforming mat
integer, intent(out) :: low
integer, intent(out) :: high

Calls

proc~~balance~~CallsGraph proc~balance balance proc~rswap RSWAP proc~balance->proc~rswap

Called by

proc~~balance~~CalledByGraph proc~balance balance proc~eigen eigen proc~eigen->proc~balance proc~colrad colrad proc~colrad->proc~eigen proc~ndmc ndmc proc~ndmc->proc~colrad proc~fida_weights_los fida_weights_los proc~fida_weights_los->proc~colrad proc~dcx dcx proc~dcx->proc~colrad proc~fida_f fida_f proc~fida_f->proc~colrad proc~pfida_f pfida_f proc~pfida_f->proc~colrad proc~read_equilibrium read_equilibrium proc~read_equilibrium->proc~colrad proc~attenuate attenuate proc~attenuate->proc~colrad proc~halo halo proc~halo->proc~colrad proc~pfida_mc pfida_mc proc~pfida_mc->proc~colrad proc~fida_weights_mc fida_weights_mc proc~fida_weights_mc->proc~colrad proc~fida_mc fida_mc proc~fida_mc->proc~colrad program~fidasim fidasim program~fidasim->proc~ndmc program~fidasim->proc~fida_weights_los program~fidasim->proc~dcx program~fidasim->proc~fida_f program~fidasim->proc~pfida_f program~fidasim->proc~read_equilibrium program~fidasim->proc~halo program~fidasim->proc~pfida_mc program~fidasim->proc~fida_weights_mc program~fidasim->proc~fida_mc proc~pnpa_mc pnpa_mc program~fidasim->proc~pnpa_mc proc~pnpa_f pnpa_f program~fidasim->proc~pnpa_f proc~npa_mc npa_mc program~fidasim->proc~npa_mc proc~npa_f npa_f program~fidasim->proc~npa_f proc~npa_weights npa_weights program~fidasim->proc~npa_weights proc~pnpa_mc->proc~attenuate proc~pnpa_f->proc~attenuate proc~npa_mc->proc~attenuate proc~npa_f->proc~attenuate proc~npa_weights->proc~attenuate

Contents

Source Code

balance

Source Code

  subroutine balance(n,  mat, scal, low, high )
    !+Balances the matrix so that the rows with zero entries
    !+off the diagonal are isolated and the remaining columns and rows
    !+are resized to have one norm close to 1.
    integer, intent(in)  :: n
        !+ Dimension of `mat`
    real(double)         :: mat(0:n-1,0:n-1)
        !+ `n`x`n` scaled matrix
    real(double)         :: scal(0:n-1)
        !+ Contains isolated eigenvalue in the positions 0-`low` and `high`-`n`-1
        !+ its other components contain the scaling factors for transforming `mat`
    integer, intent(out) :: high
    integer, intent(out) :: low


    integer, parameter   :: basis = 2
    real(double)         :: b2, r, c, f, g, s
    integer              :: m, k, i, j, iter
    scal=0.d0
    b2 = basis * basis
    m = 0
    k = n - 1
    iter=1
    do while(iter==1)
       iter = 0
       do j = k, 0, -1
          r = ZERO
          do i = 0, k
             if (i.ne.j)  r = r + DABS(mat(j,i))
          enddo
          if (r == ZERO) then
             scal(k) = j
             if (j.ne.k) then
                do i = 0, k
                   call RSWAP(mat(i,j), mat(i,k))
                enddo
                do i = m, n-1
                   call RSWAP(mat(j,i), mat(k,i))
                enddo
             endif
             k=k-1
             iter = 1
          endif
       enddo !j loop
    enddo !while iter=1
    iter=1
    do while (iter==1)
       iter = 0
       do j = m, k
          c = ZERO
          do i = m, k
             if (i.ne.j)  c = c + DABS(mat(i,j))
          enddo
          if (c == ZERO) then
             scal(m) = j
             if (j.ne.m) then
                do i = 0, k
                   call RSWAP(mat(i,j), mat(i,m))
                enddo
                do i = m, n-1
                   call RSWAP(mat(j,i), mat(m,i))
                enddo
             endif
             m = m + 1
             iter = 1
          endif
       enddo !j loop
    enddo !while iter=1
    low = m
    high = k
    do i = m, k
       scal(i) = ONE
    enddo
    iter=1
    do while (iter==1)
       iter = 0
       do i = m, k
          c=ZERO; r=ZERO
          do j = m, k
             if (j.ne.i) then
                c = c + DABS(mat(j,i))
                r = r + DABS(mat(i,j))
             endif
          enddo
          g = r / basis
          f = ONE
          s = c + r
          do while (c < g)
             f = f * basis
             c = c * b2
          enddo
          g = r * basis
          do while (c >= g)
             f = f / basis
             c = c / b2
          enddo
          if ((c + r) / f < 0.95 * s) then
             g = ONE / f
             scal(i) = scal(i) * f
             iter = 1
             do j = m, n-1
                mat(i,j) = mat(i,j) * g
             enddo
             do j = 0, k
                mat(j,i) = mat(j,i) * f
             enddo
          endif
       enddo !i loop
    enddo !while iter=1
    return
  end subroutine balance