From 873a9f9bedbbfb07d475e271923a7b86464e515f Mon Sep 17 00:00:00 2001 From: pants Date: Wed, 7 Sep 2016 14:55:30 -0400 Subject: more major refactoring --- src/net_fracture.c | 66 ++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 66 insertions(+) create mode 100644 src/net_fracture.c (limited to 'src/net_fracture.c') diff --git a/src/net_fracture.c b/src/net_fracture.c new file mode 100644 index 0000000..b7fa61d --- /dev/null +++ b/src/net_fracture.c @@ -0,0 +1,66 @@ + +#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 = alloc_break_data(net->graph->ne); + + while (true) { + double *voltages = get_voltage(net, c); + double *currents = get_current_v(net, voltages, c); + + double conductivity = get_conductivity(net, voltages, c); + + if (conductivity < 1e-12 && net->graph->boundary == TORUS_BOUND) { + free(voltages); + free(currents); + break; + } + + uint_t last_broke = get_next_broken(net, currents, cutoff); + + update_break_data(data, last_broke, fabs(conductivity * (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; +} + -- cgit v1.2.3-70-g09d2