"A heap has been corrupted" when using boost::boykov_kolmogorov_max_flow

Question

I am trying to use the Boost implementation of Boykov-Kolmogorov max-flow algorithm.

Here is how I do it:

typedef boost::adjacency_list<boost::listS, boost::vecS, boost::directedS,
boost::no_property,
boost::property<boost::edge_index_t, std::size_t> > GraphType;

typedef boost::graph_traits<GraphType>::vertex_descriptor VertexDescriptor;
typedef boost::graph_traits<GraphType>::edge_descriptor EdgeDescriptor;
typedef boost::graph_traits<GraphType>::vertices_size_type VertexIndex;
typedef boost::graph_traits<GraphType>::edges_size_type EdgeIndex;

int numberOfVertices; //The number of vertices in my graph
std::vector<std::set<int>> neighbours(numberOfVertices);

// I fill the neighbours vector with information about neighbouring vertices

int sourceId = numberOfVertices;
int sinkId = sourceId + 1;

std::size_t size = 0;
for (int i = 0; i < neighbours.size(); i++) {
    size += neighbours[i].size();
}


std::vector<float> sourceWeights(numberOfVertices);
std::vector<float> sinkWeights(numberOfVertices);

for (int i = 0; i < numberOfVertices; i++) {
    sourceWeights[i] = 1.0f;
    sinkWeights[i] = 1.0f;
}

std::vector<int> groups(numberOfVertices, 0);
std::vector<float> capacity(3 * numberOfVertices + size, 0.0f);
std::vector<EdgeDescriptor> reverseEdges(3 * numberOfVertices + size);

GraphType graph;

for (int i = 0; i < numberOfVertices; i++) {
    boost::add_vertex(graph);
}

std::vector<std::pair<int, int>> edges((size - numberOfVertices) / 2);
int countEdges = 0;
for (int i = 0; i < neighbours.size(); i++) {
    for (int c : neighbours[i]) {
        if (i < c) {
            edges[countEdges] = (std::pair<int, int>(i, c));
            countEdges++;
        }
    }
}
//I create the edges as std::pair<int,int>
//in order to add them to the graph after with the corresponding weights

for (int i = 0; i < edges.size(); i++) {
    float weight = 1.0f;

    int nextEdgeId = 2 * i;

    EdgeDescriptor edge;
    bool inserted;
    boost::tie(edge, inserted) = boost::add_edge(boost::vertex(edges[i].first, graph), boost::vertex(edges[i].second, graph), nextEdgeId, graph);
    if (!inserted)
    {
        std::cerr << "Not inserted!" << std::endl;
    }
    EdgeDescriptor reverseEdge = boost::add_edge(boost::vertex(edges[i].second, graph), boost::vertex(edges[i].first, graph), nextEdgeId + 1, graph).first;

    reverseEdges[nextEdgeId] = reverseEdge;
    reverseEdges[nextEdgeId + 1] = edge;
    capacity[nextEdgeId] = weight;
    capacity[nextEdgeId + 1] = weight;
}

VertexDescriptor sourceVertex = boost::vertex(sourceId, graph);
VertexDescriptor sinkVertex = boost::vertex(sinkId, graph);

for (int i = 0; i < numberOfVertices; i++) {
    int nextEdgeId = 2 * edges.size() + 4 * i;

    EdgeDescriptor edge;
    bool inserted;
    boost::tie(edge, inserted) = boost::add_edge(boost::vertex(i, graph), sourceVertex, nextEdgeId, graph);
    if (!inserted)
    {
        std::cerr << "Not inserted!" << std::endl;
    }
    EdgeDescriptor reverseEdge = boost::add_edge(sourceVertex, boost::vertex(i, graph), nextEdgeId + 1, graph).first;

    reverseEdges[nextEdgeId] = reverseEdge;
    reverseEdges[nextEdgeId + 1] = edge;
    capacity[nextEdgeId] = sourceWeights[i];
    capacity[nextEdgeId + 1] = sourceWeights[i];


    int nextEdgeId2 = nextEdgeId + 2;

    EdgeDescriptor edge2;
    bool inserted2;
    boost::tie(edge2, inserted2) = boost::add_edge(boost::vertex(i, graph), sinkVertex, nextEdgeId2, graph);
    if (!inserted2)
    {
        std::cerr << "Not inserted!" << std::endl;
    }
    EdgeDescriptor reverseEdge2 = boost::add_edge(sinkVertex, boost::vertex(i, graph), nextEdgeId2 + 1, graph).first;

    reverseEdges[nextEdgeId2] = reverseEdge2;
    reverseEdges[nextEdgeId2 + 1] = edge2;
    capacity[nextEdgeId2] = sinkWeights[i];
    capacity[nextEdgeId2 + 1] = sinkWeights[i];

}
std::vector<float> residual_capacity(boost::num_edges(graph), 0.0f);

//I launch the algorithm by using all the vectors I defined previously

boost::boykov_kolmogorov_max_flow(graph,
    boost::make_iterator_property_map(&capacity[0], boost::get(boost::edge_index, graph)),
    boost::make_iterator_property_map(&residual_capacity[0], boost::get(boost::edge_index, graph)),
    boost::make_iterator_property_map(&reverseEdges[0], boost::get(boost::edge_index, graph)),
    boost::make_iterator_property_map(&groups[0], boost::get(boost::vertex_index, graph)),
    boost::get(boost::vertex_index, graph),
    sourceVertex,
    sinkVertex);

I am using this code on several images (around 100). There's a lot of preprocessing on the images to extract clusters that will act as the vertices of my graph and the connectivity information between these clusters.

The code runs fine on most of the images. When I launch my program, it runs on 10 images perfectly but then, on the 11th image, for no apparent reason, it crashes with the message

Unhandled exception at 0x00007FF9B52DE6FC (ntdll.dll) in Program.exe: 0xC0000374: A heap has been corrupted (parameters: 0x00007FF9B53322B0).

It also crashes directly if I only run it on the 11th image.

The call stack is really long (thanks Boost for all those templates) so I won't post it here but the application is failing when calling boost::boykov_kolmogorov_max_flow. If I remove the line where I call this function, the program runs perfectly on all images (even the 11th one).

I am running this program in Visual Studio 2013 on a Windows 10 64-bit machine with 64Gb RAM.

So my question is, what could be the cause of this heap corruption error that only shows up on certain images (which, in this case, only means different number of vertices and neighbours information since I set the weights at 1.0 everywhere) ?