function [ParentList DiscoveryList] = Prim(EdgeListG, Weights, startV, n)

%We initialize the tree growing edge lists:
EdgeListT = []; %initially T consists of just the starting vertex startV (no edges)
%The boundary edges are all edges adjacent to startV
EdgeListGEndpts = EdgeListG(:,1:2); %indices have been removed (i.e., the third column was truncated)
[a b] = find(EdgeListGEndpts == startV); 
BdyEdges = EdgeListG(a,:);
%The boundary edges are the only edges that are initially "used"
UnusedEdges = setdiff(EdgeListG, BdyEdges, 'rows');

DiscoveryList = [startV];
ParentList = zeros(1,n);
%a zero entry in ParentList will mean the vertex has no parent "NULL"
%In the initial tree T = [startV], no vertex has a parent (and the root
%startV never will)

%Before entering into the while loop, we do the first iteration manually
%because initially EdgeListT is empty (so does not have startV)
%First Edge Selection Step (Prim)
%Need to find a lightest boundary edge
[lightest ind] = min(Weights(BdyEdges(:,3)));
e = BdyEdges(ind,:); %third entry is the index of edge w.r.t. EdgeListG
if e(1) == startV
    nonTreeEndpt = e(2);
else
   nonTreeEndpt = e(1);
end
DiscoveryList = [DiscoveryList nonTreeEndpt];
ParentList(nonTreeEndpt) = startV;
%%%Pasted from BdyEdgesUpdateInd program
%We need to delete any edge of BdyEdges that has nonTreeEndpt
%as one of its endpoints.
BdyEdgesEndpts = BdyEdges(:,1:2);
[a b] = find(BdyEdgesEndpts == nonTreeEndpt);
BdyEdges(a,:) = []; BdyEdgesEndpts(a,:) = [];
%3.  Next, we find the unused edges that are incident to nonTreeEndpt.  The
%vector a below gives their row indices in UnusedEdges, and the vector b
%gives the columnn index(1 or 2) in which nonTreeEndpt occurs.
UnusedEdgesEndpts = UnusedEdges(:,1:2);
[a b] = find(nonTreeEndpt==UnusedEdgesEndpts);
%We add these edges to BdyEdges
BdyEdges = [BdyEdges; UnusedEdges(a,:)];
%and delete them from UnusedEdges
UnusedEdges(a,:)=[];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
EdgeListT = [EdgeListT; e];

while ~isempty(BdyEdges)
    %Next Edge Selection Step (Prim)
    %Need to find a lightest boundary edge
    [lightest ind] = min(Weights(BdyEdges(:,3)));
    e = BdyEdges(ind,:); %third entry is the index of edge w.r.t. EdgeListG
    if ismember(e(1),EdgeListT(:,1:2))
        nonTreeEndpt = e(2); j = 1;
    else
        nonTreeEndpt = e(1);  j = 2;
    end
    DiscoveryList = [DiscoveryList nonTreeEndpt];
    ParentList(nonTreeEndpt) = e(j);
    [BdyEdges UnusedEdges]= BdyEdgesUpdateInd(EdgeListT, BdyEdges, e(1:2), UnusedEdges);
    EdgeListT = [EdgeListT; e];
end
%Finally we format the outputs
ParentList2 = zeros(2,n);
ParentList2(1,:) = 1:n;
ParentList2(2,:) = ParentList;
ParentList = ParentList2';