anchor.diff:      enddo
   
       num = 0
  -    100 search = 1
  -    thresh = mindeg
  -    mindeg = Neqns
  -    200 nump1 = num + 1
  -    if ( nump1>search ) search = nump1
  -    do j = search, Neqns
  -      node = Perm(j)
  -      if ( Marker(node)>=0 ) then
  -        ndeg = Deg(node)
  -        if ( ndeg<=thresh ) goto 300
  -        if ( ndeg<mindeg ) mindeg = ndeg
  -      endif
  -    enddo
  -    goto 100
  -
  -    300 search = j
  -    Nofsub = Nofsub + Deg(node)
  -    Marker(node) = 1
  -    call QMDRCH(node,Xadj,Adjncy,Deg,Marker,rchsze,Rchset,nhdsze,Nbrhd)
  -    nxnode = node
  -    do
  -      num = num + 1
  -      np = Invp(nxnode)
  -      ip = Perm(num)
  -      Perm(np) = ip
  -      Invp(ip) = np
  -      Perm(num) = nxnode
  -      Invp(nxnode) = num
  -      Deg(nxnode) = -1
  -      nxnode = Qlink(nxnode)
  -      if ( nxnode<=0 ) then
  -        if ( rchsze>0 ) then
  -          call QMDUPD(Xadj,Adjncy,rchsze,Rchset,Deg,Qsize,Qlink,Marker,Rchset(rchsze+1),Nbrhd(nhdsze+1))
  -          Marker(node) = 0
  -          do irch = 1, rchsze
  -            inode = Rchset(irch)
  -            if ( Marker(inode)>=0 ) then
  -              Marker(inode) = 0
  -              ndeg = Deg(inode)
  -              if ( ndeg<mindeg ) mindeg = ndeg
  -              if ( ndeg<=thresh ) then
  -                mindeg = thresh
  -                thresh = ndeg
  -                search = Invp(inode)
  -              endif
  +    loop1: do
  +      search = 1
  +      thresh = mindeg
  +      mindeg = Neqns
  +      loop2: do
  +        nump1 = num + 1
  +        if ( nump1>search ) search = nump1
  +        found = .false.
  +        do j = search, Neqns
  +          node = Perm(j)
  +          if ( Marker(node)>=0 ) then
  +            ndeg = Deg(node)
  +            if ( ndeg<=thresh ) then
  +              found = .true.
  +              exit
               endif
  -          enddo
  -          if ( nhdsze>0 ) call QMDOT(node,Xadj,Adjncy,Marker,rchsze,Rchset,Nbrhd)
  -        endif
  -        if ( num>=Neqns ) exit
  -        goto 200
  -      endif
  -    enddo
  +            if ( ndeg<mindeg ) mindeg = ndeg
  +          endif
  +        enddo
  +        if (.not. found) cycle loop1
  +
  +        search = j
  +        Nofsub = Nofsub + Deg(node)
  +        Marker(node) = 1
  +        call QMDRCH(node,Xadj,Adjncy,Deg,Marker,rchsze,Rchset,nhdsze,Nbrhd)
  +        nxnode = node
  +        do
  +          num = num + 1
  +          np = Invp(nxnode)
  +          ip = Perm(num)
  +          Perm(np) = ip
  +          Invp(ip) = np
  +          Perm(num) = nxnode
  +          Invp(nxnode) = num
  +          Deg(nxnode) = -1
  +          nxnode = Qlink(nxnode)
  +          if ( nxnode<=0 ) then
  +            if ( rchsze>0 ) then
  +              call QMDUPD(Xadj,Adjncy,rchsze,Rchset,Deg,Qsize,Qlink,Marker,Rchset(rchsze+1),Nbrhd(nhdsze+1))
  +              Marker(node) = 0
  +              do irch = 1, rchsze
  +                inode = Rchset(irch)
  +                if ( Marker(inode)>=0 ) then
  +                  Marker(inode) = 0
  +                  ndeg = Deg(inode)
  +                  if ( ndeg<mindeg ) mindeg = ndeg
  +                  if ( ndeg<=thresh ) then
  +                    mindeg = thresh
  +                    thresh = ndeg
  +                    search = Invp(inode)
  +                  endif
  +                endif
  +              enddo
  +              if ( nhdsze>0 ) call QMDOT(node,Xadj,Adjncy,Marker,rchsze,Rchset,Nbrhd)
  +            endif
  +            if ( num>=Neqns ) exit
  +            cycle loop2
  +          endif
  +        enddo
  +        exit
  +      enddo loop2
  +      exit
  +    enddo loop1
     end subroutine GENQMD
   
     !======================================================================!
  
positive.code:     !
      ! find zeropivot
      !
      do i = 1, Neqns
        if ( Zpiv(i)/=0 ) then
          if ( Btree(1,Invp(i))==0 ) then
            Izz = i
            goto 100
          endif
        endif
      enddo
      Izz = 0
  
      100 continue
      if ( IDBg1/=0 ) then
        write (6,"(' binary tree')")
        write (6,"(i6,'(',2I6,')')") (i,Btree(1,i),Btree(2,i),i=1,Neqns)
  
negative.code:       Parent(i) = 0
        Ancstr(i) = 0
        ip = Iperm(i)
        loop1: do k = Xadj(ip), Xadj(ip+1) - 1
          l = Invp(Adjncy(k))
          if ( l<i ) then
            do while ( Ancstr(l)/=0 )
              if ( Ancstr(l)==i ) cycle loop1
              it = Ancstr(l)
              Ancstr(l) = i
              l = it
  