GRASS GIS 7 Programmer's Manual  7.5.svn(2017)-r71817
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
centrality.c
Go to the documentation of this file.
1 /*!
2  \file lib/vector/neta/centrality.c
3 
4  \brief Network Analysis library - centrality
5 
6  Centrality measures
7 
8  (C) 2009-2010 by Daniel Bundala, and the GRASS Development Team
9 
10  This program is free software under the GNU General Public License
11  (>=v2). Read the file COPYING that comes with GRASS for details.
12 
13  \author Daniel Bundala (Google Summer of Code 2009)
14  */
15 
16 #include <stdio.h>
17 #include <stdlib.h>
18 #include <grass/gis.h>
19 #include <grass/vector.h>
20 #include <grass/glocale.h>
21 #include <grass/dgl/graph.h>
22 #include <grass/neta.h>
23 
24 /*!
25  \brief Computes degree centrality measure.
26 
27  Array degree has to be properly initialised to nnodes+1 elements
28 
29  \param graph input graph
30  \param[out] degree array of degrees
31  */
32 void NetA_degree_centrality(dglGraph_s * graph, double *degree)
33 {
34  int i;
35  double nnodes = dglGet_NodeCount(graph);
36 
37  for (i = 1; i <= nnodes; i++)
38  degree[i] =
40  dglGetNode(graph, (dglInt32_t) i)) / nnodes;
41 }
42 
43 /*!
44  \brief Computes eigenvector centrality using edge costs as weights.
45 
46  \param graph input graph
47  \param iterations number of iterations
48  \param error ?
49  \param[out] eigenvector eigen vector value
50 
51  \return 0 on success
52  \return -1 on failure
53  */
54 int NetA_eigenvector_centrality(dglGraph_s * graph, int iterations,
55  double error, double *eigenvector)
56 {
57  int i, iter, nnodes;
58  double *tmp;
59 
60  nnodes = dglGet_NodeCount(graph);
61  tmp = (double *)G_calloc(nnodes + 1, sizeof(double));
62  if (!tmp) {
63  G_fatal_error(_("Out of memory"));
64  return -1;
65  }
66 
67  error *= error;
68  for (i = 1; i <= nnodes; i++)
69  eigenvector[i] = 1;
70  for (iter = 0; iter < iterations; iter++) {
71  for (i = 1; i <= nnodes; i++)
72  tmp[i] = 0;
73  dglInt32_t *node;
76 
77  dglNode_T_Initialize(&nt, graph);
78  for (node = dglNode_T_First(&nt); node; node = dglNode_T_Next(&nt)) {
79  dglInt32_t node_id = dglNodeGet_Id(graph, node);
80  double cur_value = eigenvector[node_id];
81  dglInt32_t *edge;
82 
83  dglEdgeset_T_Initialize(&et, graph,
84  dglNodeGet_OutEdgeset(graph, node));
85  for (edge = dglEdgeset_T_First(&et); edge;
86  edge = dglEdgeset_T_Next(&et))
87  tmp[dglNodeGet_Id(graph, dglEdgeGet_Tail(graph, edge))] +=
88  cur_value * dglEdgeGet_Cost(graph, edge);
89 
91  }
92  dglNode_T_Release(&nt);
93  double cum_error = 0, max_value = tmp[1];
94 
95  for (i = 2; i <= nnodes; i++)
96  if (tmp[i] > max_value)
97  max_value = tmp[i];
98  for (i = 1; i <= nnodes; i++) {
99  tmp[i] /= max_value;
100  cum_error +=
101  (tmp[i] - eigenvector[i]) * (tmp[i] - eigenvector[i]);
102  eigenvector[i] = tmp[i];
103  }
104  if (cum_error < error)
105  break;
106 
107  }
108 
109  G_free(tmp);
110  return 0;
111 }
112 
113 /*!
114  \brief Computes betweenness and closeness centrality measure using Brandes algorithm.
115 
116  Edge costs must be nonnegative. If some edge costs are negative then
117  the behaviour of this method is undefined.
118 
119  \param graph input graph
120  \param[out] betweenness betweeness values
121  \param[out] closeness cloneness values
122 
123  \return 0 on success
124  \return -1 on failure
125  */
126 int NetA_betweenness_closeness(dglGraph_s * graph, double *betweenness,
127  double *closeness)
128 {
129  int i, j, nnodes, stack_size, count;
130  dglInt32_t *dst, *node, *stack, *cnt, *delta;
133  dglHeap_s heap;
134  struct ilist **prev;
135 
136  nnodes = dglGet_NodeCount(graph);
137 
138  dst = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
139  prev = (struct ilist **)G_calloc(nnodes + 1, sizeof(struct ilist *));
140  stack = (dglInt32_t *) G_calloc(nnodes, sizeof(dglInt32_t));
141  cnt = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
142  delta = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
143 
144  if (!dst || !prev || !stack || !cnt || !delta) {
145  G_fatal_error(_("Out of memory"));
146  return -1;
147  }
148 
149 
150  for (i = 1; i <= nnodes; i++) {
151  prev[i] = Vect_new_list();
152  if (closeness)
153  closeness[i] = 0;
154  if (betweenness)
155  betweenness[i] = 0;
156  }
157 
158  count = 0;
159  G_percent_reset();
160  dglNode_T_Initialize(&nt, graph);
161  for (node = dglNode_T_First(&nt); node; node = dglNode_T_Next(&nt)) {
162  G_percent(count++, nnodes, 1);
163  dglInt32_t s = dglNodeGet_Id(graph, node);
164  dglHeapData_u heap_data;
165  dglHeapNode_s heap_node;
166 
167  stack_size = 0;
168  for (i = 1; i <= nnodes; i++)
169  Vect_reset_list(prev[i]);
170  for (i = 1; i <= nnodes; i++) {
171  cnt[i] = 0;
172  dst[i] = -1;
173  }
174  dst[s] = 0;
175  cnt[s] = 1;
176  dglHeapInit(&heap);
177  heap_data.ul = s;
178  dglHeapInsertMin(&heap, 0, ' ', heap_data);
179  while (1) {
180  dglInt32_t v, dist;
181 
182  if (!dglHeapExtractMin(&heap, &heap_node))
183  break;
184  v = heap_node.value.ul;
185  dist = heap_node.key;
186  if (dst[v] < dist)
187  continue;
188  stack[stack_size++] = v;
189 
190  dglInt32_t *edge;
191 
192  dglEdgeset_T_Initialize(&et, graph,
193  dglNodeGet_OutEdgeset(graph,
194  dglGetNode(graph,
195  v)));
196  for (edge = dglEdgeset_T_First(&et); edge;
197  edge = dglEdgeset_T_Next(&et)) {
198  dglInt32_t *to = dglEdgeGet_Tail(graph, edge);
199  dglInt32_t to_id = dglNodeGet_Id(graph, to);
200  dglInt32_t d = dglEdgeGet_Cost(graph, edge);
201 
202  if (dst[to_id] == -1 || dst[to_id] > dist + d) {
203  dst[to_id] = dist + d;
204  Vect_reset_list(prev[to_id]);
205  heap_data.ul = to_id;
206  dglHeapInsertMin(&heap, dist + d, ' ', heap_data);
207  }
208  if (dst[to_id] == dist + d) {
209  cnt[to_id] += cnt[v];
210  Vect_list_append(prev[to_id], v);
211  }
212  }
213 
215  }
216  dglHeapFree(&heap, NULL);
217  for (i = 1; i <= nnodes; i++)
218  delta[i] = 0;
219  for (i = stack_size - 1; i >= 0; i--) {
220  dglInt32_t w = stack[i];
221 
222  if (closeness)
223  closeness[s] += dst[w];
224 
225  for (j = 0; j < prev[w]->n_values; j++) {
226  dglInt32_t v = prev[w]->value[j];
227 
228  delta[v] += (cnt[v] / (double)cnt[w]) * (1.0 + delta[w]);
229  }
230  if (w != s && betweenness)
231  betweenness[w] += delta[w];
232 
233  }
234  if (closeness)
235  closeness[s] /= (double)stack_size;
236  }
237  dglNode_T_Release(&nt);
238 
239  for (i = 1; i <= nnodes; i++)
240  Vect_destroy_list(prev[i]);
241  G_free(delta);
242  G_free(cnt);
243  G_free(stack);
244  G_free(prev);
245  G_free(dst);
246 
247  return 0;
248 };
void G_free(void *buf)
Free allocated memory.
Definition: gis/alloc.c:149
real i
Definition: la.h:56
dglInt32_t * dglGetNode(dglGraph_s *pGraph, dglInt32_t nNodeId)
long key
Definition: heap.h:38
void Vect_destroy_list(struct ilist *list)
Frees all memory associated with a struct ilist, including the struct itself.
dglHeapData_u value
Definition: heap.h:39
struct ilist * Vect_new_list(void)
Creates and initializes a struct ilist.
dglInt32_t dglNodeGet_Id(dglGraph_s *pGraph, dglInt32_t *pnNode)
void dglNode_T_Release(dglNodeTraverser_s *pT)
dglInt32_t dglEdgeGet_Cost(dglGraph_s *pGraph, dglInt32_t *pnEdge)
void dglHeapFree(dglHeap_s *pheap, dglHeapCancelItem_fn pfnCancelItem)
Definition: heap.c:37
int n_values
Number of values in the list.
Definition: gis.h:659
int count
Definition: heap.h:44
char * dst
Definition: lz4.h:354
int dglGet_NodeCount(dglGraph_s *pgraph)
void dglHeapInit(dglHeap_s *pheap)
Definition: heap.c:29
dglInt32_t * dglEdgeGet_Tail(dglGraph_s *pGraph, dglInt32_t *pnEdge)
#define NULL
Definition: ccmath.h:32
void G_fatal_error(const char *msg,...)
Print a fatal error message to stderr.
Definition: gis/error.c:159
dglInt32_t * dglNode_T_First(dglNodeTraverser_s *pT)
dglInt32_t * dglEdgeset_T_First(dglEdgesetTraverser_s *pT)
int Vect_reset_list(struct ilist *list)
Reset ilist structure.
long dglInt32_t
Definition: type.h:37
int Vect_list_append(struct ilist *list, int val)
Append new item to the end of list if not yet present.
void G_percent_reset(void)
Reset G_percent() to 0%; do not add newline.
Definition: percent.c:119
void dglEdgeset_T_Release(dglEdgesetTraverser_s *pT)
void G_percent(long n, long d, int s)
Print percent complete messages.
Definition: percent.c:62
int NetA_eigenvector_centrality(dglGraph_s *graph, int iterations, double error, double *eigenvector)
Computes eigenvector centrality using edge costs as weights.
Definition: centrality.c:54
void NetA_degree_centrality(dglGraph_s *graph, double *degree)
Computes degree centrality measure.
Definition: centrality.c:32
int dglNodeGet_OutDegree(dglGraph_s *pGraph, dglInt32_t *pnNode)
#define _(str)
Definition: glocale.h:13
List of integers.
Definition: gis.h:650
int * value
Array of values.
Definition: gis.h:655
int dglHeapInsertMin(dglHeap_s *pheap, long key, unsigned char flags, dglHeapData_u value)
Definition: heap.c:52
int dglHeapExtractMin(dglHeap_s *pheap, dglHeapNode_s *pnoderet)
Definition: heap.c:79
int dglNode_T_Initialize(dglNodeTraverser_s *pT, dglGraph_s *pGraph)
dglInt32_t * dglEdgeset_T_Next(dglEdgesetTraverser_s *pT)
dglInt32_t * dglNodeGet_OutEdgeset(dglGraph_s *pGraph, dglInt32_t *pnNode)
int dglEdgeset_T_Initialize(dglEdgesetTraverser_s *pT, dglGraph_s *pGraph, dglInt32_t *pnEdgeset)
dglInt32_t * dglNode_T_Next(dglNodeTraverser_s *pT)
int NetA_betweenness_closeness(dglGraph_s *graph, double *betweenness, double *closeness)
Computes betweenness and closeness centrality measure using Brandes algorithm.
Definition: centrality.c:126
unsigned long ul
Definition: heap.h:31