GRASS Programmer's Manual  6.4.2(2012)
centrality.c
Go to the documentation of this file.
00001 
00016 #include <stdio.h>
00017 #include <stdlib.h>
00018 #include <grass/gis.h>
00019 #include <grass/Vect.h>
00020 #include <grass/glocale.h>
00021 #include <grass/dgl/graph.h>
00022 #include <grass/neta.h>
00023 
00032 void NetA_degree_centrality(dglGraph_s * graph, double *degree)
00033 {
00034     int i;
00035     double nnodes = dglGet_NodeCount(graph);
00036 
00037     for (i = 1; i <= nnodes; i++)
00038         degree[i] =
00039             dglNodeGet_OutDegree(graph,
00040                                  dglGetNode(graph, (dglInt32_t) i)) / nnodes;
00041 }
00042 
00054 int NetA_eigenvector_centrality(dglGraph_s * graph, int iterations,
00055                                 double error, double *eigenvector)
00056 {
00057     int i, iter, nnodes;
00058     double *tmp;
00059 
00060     nnodes = dglGet_NodeCount(graph);
00061     tmp = (double *)G_calloc(nnodes + 1, sizeof(double));
00062     if (!tmp) {
00063         G_fatal_error(_("Out of memory"));
00064         return -1;
00065     }
00066 
00067     error *= error;
00068     for (i = 1; i <= nnodes; i++)
00069         eigenvector[i] = 1;
00070     for (iter = 0; iter < iterations; iter++) {
00071         for (i = 1; i <= nnodes; i++)
00072             tmp[i] = 0;
00073         dglInt32_t *node;
00074         dglNodeTraverser_s nt;
00075         dglEdgesetTraverser_s et;
00076 
00077         dglNode_T_Initialize(&nt, graph);
00078         for (node = dglNode_T_First(&nt); node; node = dglNode_T_Next(&nt)) {
00079             dglInt32_t node_id = dglNodeGet_Id(graph, node);
00080             double cur_value = eigenvector[node_id];
00081             dglInt32_t *edge;
00082 
00083             dglEdgeset_T_Initialize(&et, graph,
00084                                     dglNodeGet_OutEdgeset(graph, node));
00085             for (edge = dglEdgeset_T_First(&et); edge;
00086                  edge = dglEdgeset_T_Next(&et))
00087                 tmp[dglNodeGet_Id(graph, dglEdgeGet_Tail(graph, edge))] +=
00088                     cur_value * dglEdgeGet_Cost(graph, edge);
00089 
00090             dglEdgeset_T_Release(&et);
00091         }
00092         dglNode_T_Release(&nt);
00093         double cum_error = 0, max_value = tmp[1];
00094 
00095         for (i = 2; i <= nnodes; i++)
00096             if (tmp[i] > max_value)
00097                 max_value = tmp[i];
00098         for (i = 1; i <= nnodes; i++) {
00099             tmp[i] /= max_value;
00100             cum_error +=
00101                 (tmp[i] - eigenvector[i]) * (tmp[i] - eigenvector[i]);
00102             eigenvector[i] = tmp[i];
00103         }
00104         if (cum_error < error)
00105             break;
00106 
00107     }
00108 
00109     G_free(tmp);
00110     return 0;
00111 }
00112 
00126 int NetA_betweenness_closeness(dglGraph_s * graph, double *betweenness,
00127                                double *closeness)
00128 {
00129     int i, j, nnodes, stack_size, count;
00130     dglInt32_t *dst, *node, *stack, *cnt, *delta;
00131     dglNodeTraverser_s nt;
00132     dglEdgesetTraverser_s et;
00133     dglHeap_s heap;
00134     struct ilist **prev;
00135 
00136     nnodes = dglGet_NodeCount(graph);
00137 
00138     dst = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
00139     prev = (struct ilist **)G_calloc(nnodes + 1, sizeof(struct ilist *));
00140     stack = (dglInt32_t *) G_calloc(nnodes, sizeof(dglInt32_t));
00141     cnt = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
00142     delta = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
00143 
00144     if (!dst || !prev || !stack || !cnt || !delta) {
00145         G_fatal_error(_("Out of memory"));
00146         return -1;
00147     }
00148 
00149 
00150     for (i = 1; i <= nnodes; i++) {
00151         prev[i] = Vect_new_list();
00152         if (closeness)
00153             closeness[i] = 0;
00154         if (betweenness)
00155             betweenness[i] = 0;
00156     }
00157 
00158     count = 0;
00159     G_percent_reset();
00160     dglNode_T_Initialize(&nt, graph);
00161     for (node = dglNode_T_First(&nt); node; node = dglNode_T_Next(&nt)) {
00162         G_percent(count++, nnodes, 1);
00163         dglInt32_t s = dglNodeGet_Id(graph, node);
00164         dglHeapData_u heap_data;
00165         dglHeapNode_s heap_node;
00166 
00167         stack_size = 0;
00168         for (i = 1; i <= nnodes; i++)
00169             Vect_reset_list(prev[i]);
00170         for (i = 1; i <= nnodes; i++) {
00171             cnt[i] = 0;
00172             dst[i] = -1;
00173         }
00174         dst[s] = 0;
00175         cnt[s] = 1;
00176         dglHeapInit(&heap);
00177         heap_data.ul = s;
00178         dglHeapInsertMin(&heap, 0, ' ', heap_data);
00179         while (1) {
00180             dglInt32_t v, dist;
00181 
00182             if (!dglHeapExtractMin(&heap, &heap_node))
00183                 break;
00184             v = heap_node.value.ul;
00185             dist = heap_node.key;
00186             if (dst[v] < dist)
00187                 continue;
00188             stack[stack_size++] = v;
00189 
00190             dglInt32_t *edge;
00191 
00192             dglEdgeset_T_Initialize(&et, graph,
00193                                     dglNodeGet_OutEdgeset(graph,
00194                                                           dglGetNode(graph,
00195                                                                      v)));
00196             for (edge = dglEdgeset_T_First(&et); edge;
00197                  edge = dglEdgeset_T_Next(&et)) {
00198                 dglInt32_t *to = dglEdgeGet_Tail(graph, edge);
00199                 dglInt32_t to_id = dglNodeGet_Id(graph, to);
00200                 dglInt32_t d = dglEdgeGet_Cost(graph, edge);
00201 
00202                 if (dst[to_id] == -1 || dst[to_id] > dist + d) {
00203                     dst[to_id] = dist + d;
00204                     Vect_reset_list(prev[to_id]);
00205                     heap_data.ul = to_id;
00206                     dglHeapInsertMin(&heap, dist + d, ' ', heap_data);
00207                 }
00208                 if (dst[to_id] == dist + d) {
00209                     cnt[to_id] += cnt[v];
00210                     Vect_list_append(prev[to_id], v);
00211                 }
00212             }
00213 
00214             dglEdgeset_T_Release(&et);
00215         }
00216         dglHeapFree(&heap, NULL);
00217         for (i = 1; i <= nnodes; i++)
00218             delta[i] = 0;
00219         for (i = stack_size - 1; i >= 0; i--) {
00220             dglInt32_t w = stack[i];
00221 
00222             if (closeness)
00223                 closeness[s] += dst[w];
00224 
00225             for (j = 0; j < prev[w]->n_values; j++) {
00226                 dglInt32_t v = prev[w]->value[j];
00227 
00228                 delta[v] += (cnt[v] / (double)cnt[w]) * (1.0 + delta[w]);
00229             }
00230             if (w != s && betweenness)
00231                 betweenness[w] += delta[w];
00232 
00233         }
00234         if (closeness)
00235             closeness[s] /= (double)stack_size;
00236     }
00237     dglNode_T_Release(&nt);
00238 
00239     for (i = 1; i <= nnodes; i++)
00240         Vect_destroy_list(prev[i]);
00241     G_free(delta);
00242     G_free(cnt);
00243     G_free(stack);
00244     G_free(prev);
00245     G_free(dst);
00246 
00247     return 0;
00248 };
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines