allpairs.c

Go to the documentation of this file.

#include <stdio.h>
#include <stdlib.h>
#include <grass/gis.h>
#include <grass/Vect.h>
#include <grass/glocale.h>
#include <grass/dgl/graph.h>

int NetA_allpairs(dglGraph_s * graph, dglInt32_t ** dist)
{
    int nnodes, i, j, k, indices;
    dglEdgesetTraverser_s et;
    dglNodeTraverser_s nt;
    dglInt32_t *node;

    nnodes = dglGet_NodeCount(graph);
    dglInt32_t *node_indices;

    node_indices = (dglInt32_t *) G_calloc(nnodes, sizeof(dglInt32_t));
    if (!node_indices) {
        G_fatal_error(_("Out of memory"));
        return -1;
    }
    G_message(_("Computing all pairs shortest paths..."));
    G_percent_reset();
    for (i = 0; i <= nnodes; i++)
        for (j = 0; j <= nnodes; j++)
            dist[i][j] = -1;
    dglNode_T_Initialize(&nt, graph);
    indices = 0;
    for (node = dglNode_T_First(&nt); node; node = dglNode_T_Next(&nt)) {
        dglInt32_t node_id = dglNodeGet_Id(graph, node);

        node_indices[indices++] = node_id;
        dglInt32_t *edge;

        dglEdgeset_T_Initialize(&et, graph,
                                dglNodeGet_OutEdgeset(graph, node));
        for (edge = dglEdgeset_T_First(&et); edge;
             edge = dglEdgeset_T_Next(&et))
            if (dglEdgeGet_Id(graph, edge) < 0) /*ignore backward edges */
                dist[node_id][dglNodeGet_Id
                              (graph, dglEdgeGet_Tail(graph, edge))] =
                    dglEdgeGet_Cost(graph, edge);
        dglEdgeset_T_Release(&et);
    }
    dglNode_T_Release(&nt);
    for (k = 0; k < indices; k++) {
        dglInt32_t k_index = node_indices[k];

        G_percent(k + 1, indices, 1);
        for (i = 0; i < indices; i++) {
            dglInt32_t i_index = node_indices[i];

            if (dist[i_index][k_index] == -1)
                continue;       /*no reason to proceed along infinite path */
            for (j = 0; j < indices; j++) {
                dglInt32_t j_index = node_indices[j];

                if (dist[k_index][j_index] != -1 &&
                    (dist[i_index][k_index] + dist[k_index][j_index] <
                     dist[i_index][j_index] ||
                     dist[i_index][j_index] == -1)) {
                    dist[i_index][j_index] =
                        dist[i_index][k_index] + dist[k_index][j_index];
                }
            }
        }
    }

    G_free(node_indices);
    return 0;
}