1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
|
#include "fracture.h"
uint_t get_next_broken(net_t *net, double *currents, double cutoff) {
uint_t max_pos = UINT_MAX;
double max_val = 0;
for (uint_t i = 0; i < net->graph->ne; i++) {
double current = fabs(currents[i]);
bool broken = net->fuses[i];
if (!broken && current > cutoff) {
double val = current / net->thres[i];
if (val > max_val) {
max_val = val;
max_pos = i;
}
}
}
if (max_pos == UINT_MAX) {
printf("GET_NEXT_BROKEN: currents is zero or NaN, no max_val found\n");
exit(EXIT_FAILURE);
}
return max_pos;
}
data_t *net_fracture(net_t *net, cholmod_common *c, double cutoff) {
data_t *data = data_create(net->graph->ne);
while (true) {
double *voltages = net_voltages(net, c);
double *currents = net_currents(net, voltages, c);
double conductivity = net_conductivity(net, voltages);
if (conductivity < 1e-12 && net->graph->boundary == TORUS_BOUND) {
free(voltages);
free(currents);
break;
}
uint_t last_broke = get_next_broken(net, currents, cutoff);
double sim_current;
if (net->voltage_bound) {
sim_current = conductivity;
} else {
sim_current = 1;
}
data_update(data, last_broke, fabs(sim_current * (net->thres)[last_broke] / currents[last_broke]), conductivity);
free(voltages);
free(currents);
break_edge(net, last_broke, c);
if (net->num_components > 1 && net->graph->boundary == TORUS_BOUND) {
break;
}
if (net->marks[net->graph->nv] != net->marks[net->graph->nv + 1] && net->graph->boundary != TORUS_BOUND) {
break;
}
}
return data;
}
|