generate_adjacency Module Subroutine

  module subroutine generate_adjacency(this, index_list)
    !! Generate the adjacency matrix of the graph.
    implicit none

    ! Arguments
    class(graph_type), intent(inout) :: this
    !! Parent. Instance of the graph structure.
    integer, dimension(:,:), intent(in), optional :: index_list
    !! List of indices to be used for the adjacency matrix.

    integer :: i, j, k, i1, i2, num_edges
    !! Loop indices.
    integer, dimension(3,this%num_vertices+2*this%num_edges) :: coo
    !! sparse adjacency in coordinate format (coo)
    integer, dimension(3) :: temp
    !! temporary array

    if( .not. this%is_sparse )then  ! graph is not sparse; allocate and fill array adjacency.
       if(allocated(this%adjacency)) deallocate(this%adjacency)
       allocate(this%adjacency(this%num_vertices, this%num_vertices))
       this%adjacency = 0
       do k = 1, this%num_edges
          i = this%edge(k)%index(1)
          j = this%edge(k)%index(2)
          if(i.eq.0.or.j.eq.0) cycle  ! skip edges with zero index
          if(this%directed.and.j.lt.0)then
             this%adjacency(i,abs(j)) = k
          else
             this%adjacency(i,abs(j)) = k
             this%adjacency(abs(j),i) = k
          end if
       end do
    else
       if(allocated(this%adj_ia)) deallocate(this%adj_ia)
       if(allocated(this%adj_ja)) deallocate(this%adj_ja)
       allocate(this%adj_ia(this%num_vertices+1))
       if(this%directed)then
          num_edges = this%num_edges
       else
          num_edges = 2*this%num_edges
          if(this%has_self_loops) num_edges = num_edges - this%num_vertices
       end if
       allocate(this%adj_ja(2, num_edges))
       ! Step 1: edgelist to coo array.
       !  do i = 1, this%num_vertices  ! vertex self-loops
       !     coo(1,i) = i
       !     coo(2,i) = i
       !     coo(3,i) = 0
       !  end do
       if(present(index_list))then
          do i = 1, size(index_list, dim=2)
             coo(1,i) = index_list(1,i)
             coo(2,i) = index_list(2,i)
             coo(3,i) = i
          end do
          j = this%num_edges
          if(.not.this%directed)then
             do i = 1, size(index_list, dim=2)
                if(index_list(1,i).eq.index_list(2,i)) cycle ! skip self-loops
                j = j + 1
                coo(1,j) = index_list(2,i)
                coo(2,j) = index_list(1,i)
                coo(3,j) = i
             end do
          end if
       else
          do i = 1, this%num_edges     ! first pass over edges (1,2)
             coo(1,i) = this%edge(i)%index(1)
             coo(2,i) = this%edge(i)%index(2)
             coo(3,i) = i
          end do
          if(.not.this%directed)then
             j = this%num_edges
             do i = 1, this%num_edges     ! second pass over edges (2,1)
                if(this%edge(i)%index(1).eq.this%edge(i)%index(2)) cycle ! skip self-loops
                j = j + 1
                coo(1,j) = this%edge(i)%index(2)
                coo(2,j) = this%edge(i)%index(1)
                coo(3,j) = i
             end do
          end if
       end if
       ! Step 2: sort coo array.
       do i2 = num_edges, 2, -1
          do i1 = 1, i2 - 1
             if( ( coo(1,i1) .gt. ( coo(1,i2) ) ) .or. &
                  ( &
                       ( coo(1,i1) .eq. coo(1,i2) ) .and. &
                       ( coo(2,i1) .gt. coo(2,i2)  ) &
                  ) &
             )then
                temp(1:3) = coo(1:3,i1)
                coo(1:3,i1) = coo(1:3,i2)
                coo(1:3,i2) = temp(1:3)
             end if
          end do
       end do
       ! Step 3: sorted coo array to adj_ia and adj_ja arrays.
       this%adj_ia(1) = 1
       do i = 1, num_edges
          this%adj_ja(1,i) = coo(2,i)
          this%adj_ja(2,i) = coo(3,i)
          this%adj_ia(coo(1,i)+1) = i + 1
       end do
    end if

  end subroutine generate_adjacency