anchor.diff:                     MPI_COMM_WORLD,ierr)                                 !
       call MPI_BCAST(A%val,A%nnz,MPI_fkind,D_MASTER,&                     !
                      MPI_COMM_WORLD,ierr)                                 !
  -    call SpMtx_arrange_clrorder(A,numprocs,M%eptnmap,partorder, &       !
  -                                partstarts,sort=.true.)                 !
  +    write(stream,*)'A orig:'
  +    call SpMtx_printRaw(A)
  +    call SpMtx_arrange(A,arrange_type=D_SpMtx_ARRNG_ROWS,sort=.true.)   !
  +    !========= count color elements ============
  +    if (A%arrange_type==D_SpMtx_ARRNG_ROWS) then
  +      n=A%nrows
  +    elseif (A%arrange_type==D_SpMtx_ARRNG_COLS) then
  +      n=A%ncols
  +    endif
  +
  +    allocate(ccount(numprocs))
  +    ccount=0
  +    do i=1,n
  +      ccount(M%eptnmap(i))=ccount(M%eptnmap(i))+1
  +    enddo
  +    allocate(clrstarts(numprocs+1))
  +    clrstarts(1)=1
  +    do i=1,numprocs
  +      clrstarts(i+1)=clrstarts(i)+ccount(i)
  +    end do
  +    allocate(clrorder(n))
  +    ccount(1:numprocs)=clrstarts(1:numprocs)
  +    do i=1,n
  +      clrorder(ccount(M%eptnmap(i)))=i
  +      ccount(M%eptnmap(i))=ccount(M%eptnmap(i))+1
  +    enddo
       do i=1,numprocs                                                     !
         write(stream,*)'partition ',i,' is in:', &                        !
  -         partorder(partstarts(i):partstarts(i+1)-1)                     !
  +         clrorder(clrstarts(i):clrstarts(i+1)-1)                     !
       enddo                                                               !
  +    deallocate(ccount)
       !-------------------------------------------------------------------+
  -    call SpMtx_keep_subd_wol(myrank+1,max(sctls%overlap,sctls%smoothers),&
  -                             A,M,partorder,partstarts) 
  -    !call SpMtx_printRaw(A)
  +    write(stream,*)'A after arrange:'
  +    call SpMtx_printRaw(A)
  +    call SpMtx_build_interfghost(myrank+1,max(sctls%overlap,sctls%smoothers),&
  +                             A,A_interf,A_ghost,M,clrorder,clrstarts) 
  +    write(stream,*)'A after:'
  +    call SpMtx_printRaw(A)
  +    write(stream,*)'A%M_bound=',(i,':',A%M_bound(i),' ',i=1,A%nrows+1)
  +    write(stream,*)'A interf:'
  +    call SpMtx_printRaw(A_interf)
  +    write(stream,*)'A ghost:'
  +    call SpMtx_printRaw(A_ghost)
       call MPI_BARRIER(MPI_COMM_WORLD,ierr)
       call DOUG_abort('testing nodal graph partitioning',0)
     end subroutine SpMtx_DistributeAssembled              
   
     ! Take away from matrix unneeded elements...
     ! (the matrix should be arranged into row format with SpMtx_arrange_clrorder)
  -  subroutine SpMtx_keep_subd_wol(clr,ol,A,M,clrorder,clrstarts)
  -    integer,intent(in)                 :: clr      !color #
  +  subroutine SpMtx_build_interfghost(clr,ol,A,A_interf,A_ghost,M,clrorder,clrstarts)
  +    integer,intent(in)                 :: clr      !the color # we are keeping
       integer,intent(in)                 :: ol       !overlap size
       Type(SpMtx), intent(in out)        :: A        !Initial matrix
  +    Type(SpMtx), intent(in out)        :: A_interf !new mat.betw. A and A_ghost
  +    Type(SpMtx), intent(in out)        :: A_ghost  !matrix on ghost nodes
       type(Mesh)                         :: M        !Mesh object
       integer,dimension(:),pointer       :: clrorder
        !order of matrix rows (columns) so that color i is found in rows (columns):
                      !             clrorder(clrstarts(i):clrstarts(i+1)-1)
  -    integer,dimension(:),pointer       :: clrstarts  !allocated here
  +    integer,dimension(:),pointer       :: clrstarts  !(allocated earlier)
       !local:
  -    integer :: i,j,clrnode,clrneigh,nfront
  +    integer :: i,j,k,clrnode,clrneigh,nfront,layer,lastlayer,node,nnz
  +    integer :: maxclrnode
       integer,dimension(:),pointer :: helper,front
       integer,dimension(:),pointer :: onfront
  +    integer,dimension(:),pointer :: frontstart,frontend
  +    integer :: ghostsz,interfsz
  +    integer,dimension(:),pointer :: itmp,jtmp,btmp
  +    float(kind=rk), dimension(:), pointer :: rtmp
  +
       
       !write(stream,*)'BIT_SIZE(1)',BIT_SIZE(helper(2))
       allocate(helper(numprocs))
       helper=0
  -    allocate(front(A%nrows)) ! for keeping track on the front
  +    allocate(front(A%nrows)) ! for keeping track of the front
       allocate(onfront(A%nrows)) ! for keeping track on the front
  +    onfront=0
  +    allocate(frontstart(0:ol))
  +    allocate(frontend(0:ol))
       nfront=0
  +    frontstart(0)=1
  +    interfsz=0
       !first, count how many neighbours the clr has
       if (A%arrange_type==D_SpMtx_ARRNG_ROWS) then
         do i=clrstarts(clr),clrstarts(clr+1)-1
           clrnode=clrorder(i)
           do j=A%M_bound(clrnode),A%M_bound(clrnode+1)-1
  -          clrneigh=M%eptnmap(A%indj(j))
  +          node=A%indj(j)
  +          clrneigh=M%eptnmap(node)
             if (clrneigh/=clr) then
               helper(clrneigh)=1
  -            if (onfront(A%indj(j))==0) then
  +            interfsz=interfsz+1
  +            if (onfront(clrnode)==0) then
                 nfront=nfront+1
  -              front(nfront)=A%indj(j)
  -              onfront(A%indj(j))=1
  +              front(nfront)=clrnode
  +              onfront(clrnode)=-1
               endif
             endif
           enddo
  
positive.code:     endif
    end subroutine SpMtx_aggregate
  
    subroutine SpMtx_DistributeAssembled(A,M)
      use Graph_class
      use Mesh_class
      implicit none
  
      type(SpMtx),intent(inout) :: A
      type(Mesh)                :: M
      integer :: i,ierr
      integer, dimension(:), pointer :: xadj
      integer, dimension(:), pointer :: adjncy
      integer                        :: nedges
      type(Graph) :: G
      integer, dimension(6) :: part_opts = (/0,0,0,0,0,0/)
      integer,dimension(:),pointer       :: partorder,partstarts
      integer,dimension(4)           :: buf
      
      if (ismaster()) then ! Here master simply splits the matrix into pieces
  
negative.code:     else
        call DOUG_abort('SpMtx_keep_subd_wol: Matrix arrangment not done!',19)
      endif
      frontend(0)=nfront
      !write(stream,*)'I have ',sum(helper),' neighbours!'
      M%nnghbrs=sum(helper)
      allocate(M%nghbrs(M%nnghbrs))
  