From 03de79b8c5ebcc206e3450dfbc701211d9c254b0 Mon Sep 17 00:00:00 2001 From: pants Date: Fri, 9 Sep 2016 14:33:56 -0400 Subject: more refactoring --- src/net.c | 130 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 130 insertions(+) create mode 100644 src/net.c (limited to 'src/net.c') diff --git a/src/net.c b/src/net.c new file mode 100644 index 0000000..9f6965b --- /dev/null +++ b/src/net.c @@ -0,0 +1,130 @@ + +#include "fracture.h" + +double *get_thres(uint_t ne, double beta) { + double *thres = (double *)malloc(ne * sizeof(double)); + assert(thres != NULL); + + gsl_rng *r = gsl_rng_alloc(GSL_RAND_GEN); + gsl_rng_set(r, rand_seed()); + + for (uint_t i = 0; i < ne; i++) { + thres[i] = rand_dist_pow(r, beta); + } + + gsl_rng_free(r); + + return thres; +} + +void net_notch(net_t *net, double notch_len, cholmod_common *c) { + for (uint_t i = 0; i < net->graph->ne; i++) { + uint_t v1, v2; + double v1x, v1y, v2x, v2y, dy; + bool crosses_center, not_wrapping, correct_length; + + v1 = net->graph->ev[2 * i]; + v2 = net->graph->ev[2 * i + 1]; + + v1x = net->graph->vx[2 * v1]; + v1y = net->graph->vx[2 * v1 + 1]; + v2x = net->graph->vx[2 * v2]; + v2y = net->graph->vx[2 * v2 + 1]; + + dy = v1y - v2y; + + crosses_center = (v1y >= 0.5 && v2y <= 0.5) || (v1y <= 0.5 && v2y >= 0.5); + not_wrapping = fabs(dy) < 0.5; + //correct_length = v1x + dx / dy * (v1y - 0.5) <= notch_len; + correct_length = v1x < notch_len && v2x < notch_len; + + if (crosses_center && not_wrapping && correct_length) { + break_edge(net, i, c); + } + } +} + +net_t *net_create(const graph_t *g, double inf, double beta, double notch_len, bool vb, cholmod_common *c) { + net_t *net = (net_t *)calloc(1, sizeof(net_t)); + assert(net != NULL); + + net->graph = g; + net->fuses = (bool *)calloc(g->ne, sizeof(bool)); + assert(net->fuses != NULL); + net->thres = get_thres(g->ne, beta); + net->inf = inf; + + net->voltage_bound = vb; + net->boundary_cond = bound_set(g, vb, notch_len, c); + + if (g->boundary != TORUS_BOUND) net->adjacency = gen_adjacency(net, false, false, 0, c); + else net->adjacency = gen_adjacency(net, true, false, 0, c); + + net->marks = (uint_t *)malloc((net->graph->break_dim) * sizeof(uint_t)); + assert(net->marks != NULL); + + net->dual_marks = (uint_t *)malloc((net->graph->dnv) * sizeof(uint_t)); + assert(net->dual_marks != NULL); + + for (uint_t i = 0; i < (net->graph->break_dim); i++) { + net->marks[i] = 1; + } + for (uint_t i = 0; i < (net->graph->dnv); i++) { + net->dual_marks[i] = i+1; + } + net->num_components = 1; + + net_notch(net, notch_len, c); + + { + cholmod_sparse *laplacian = gen_laplacian(net, c, true); + net->factor = CHOL_F(analyze)(laplacian, c); + CHOL_F(factorize)(laplacian, net->factor, c); + CHOL_F(free_sparse)(&laplacian, c); + } + + return net; +} + +net_t *net_copy(const net_t *net, cholmod_common *c) { + net_t *net_copy = (net_t *)calloc(1, sizeof(net_t)); + assert(net_copy != NULL); + memcpy(net_copy, net, sizeof(net_t)); + + size_t fuses_size = (net->graph)->ne * sizeof(bool); + net_copy->fuses = (bool *)malloc(fuses_size); + assert(net_copy->fuses != NULL); + memcpy(net_copy->fuses, net->fuses, fuses_size); + + size_t thres_size = (net->graph)->ne * sizeof(double); + net_copy->thres = (double *)malloc(thres_size); + assert(net_copy->thres != NULL); + memcpy(net_copy->thres, net->thres, thres_size); + + size_t marks_size = (net->graph->break_dim) * sizeof(uint_t); + net_copy->marks = (uint_t *)malloc(marks_size); + assert(net_copy->marks != NULL); + memcpy(net_copy->marks, net->marks, marks_size); + + size_t dual_marks_size = (net->graph->dnv) * sizeof(uint_t); + net_copy->dual_marks = (uint_t *)malloc(dual_marks_size); + assert(net_copy->dual_marks != NULL); + memcpy(net_copy->dual_marks, net->dual_marks, dual_marks_size); + + net_copy->adjacency = CHOL_F(copy_sparse)(net->adjacency, c); + net_copy->boundary_cond = CHOL_F(copy_dense)(net->boundary_cond, c); + net_copy->factor = CHOL_F(copy_factor)(net->factor, c); + + return net_copy; +} + +void net_free(net_t *net, cholmod_common *c) { + free(net->fuses); + free(net->thres); + CHOL_F(free_dense)(&(net->boundary_cond), c); + CHOL_F(free_sparse)(&(net->adjacency), c); + CHOL_F(free_factor)(&(net->factor), c); + free(net->marks); + free(net->dual_marks); + free(net); +} -- cgit v1.2.3-70-g09d2