#include "fracture.h" double *gen_fuse_thres(unsigned int num_edges, double *edge_coords, double beta, double (*beta_scale)(double, double, double)) { unsigned int size = num_edges; assert(beta > 0); double *fuse_thres = (double *)malloc(size * sizeof(double)); assert(fuse_thres != NULL); gsl_set_error_handler_off(); gsl_rng *r = gsl_rng_alloc(gsl_rng_mt19937); FILE *rf = fopen("/dev/urandom", "r"); unsigned long int seed; fread(&seed, sizeof(unsigned long int), 1, rf); fclose(rf); gsl_rng_set(r, seed); for (unsigned int i = 0; i < size; i++) { double x = edge_coords[2 * i]; double y = edge_coords[2 * i + 1]; while ((fuse_thres[i] = gsl_sf_exp(gsl_sf_log(gsl_ran_flat(r, 0, 1)) / beta_scale(beta, x, y))) == 0.0) ; } gsl_rng_free(r); return fuse_thres; } void gen_voro_crack(net_t *instance, double crack_len, cholmod_common *c) { for (uint_t i = 0; i < instance->graph->ne; i++) { uint_t v1, v2; double v1x, v1y, v2x, v2y, dx, dy; v1 = instance->graph->ev[2 * i]; v2 = instance->graph->ev[2 * i + 1]; v1x = instance->graph->vert_coords[2 * v1]; v1y = instance->graph->vert_coords[2 * v1 + 1]; v2x = instance->graph->vert_coords[2 * v2]; v2y = instance->graph->vert_coords[2 * v2 + 1]; dx = v1x - v2x; dy = v1y - v2y; if (((v1y > 0.5 && v2y < 0.5) || (v1y < 0.5 && v2y > 0.5)) && fabs(dy) < 0.5) { if (v1x + dx / dy * (v1y - 0.5) <= crack_len) { break_edge(instance, i, c); } } } } bool gen_crack(net_t *instance, double crack_len, double crack_width, cholmod_common *c) { assert(instance != NULL); bool *fuses = instance->fuses; assert(fuses != NULL); const graph_t *network = instance->graph; unsigned int num_edges = network->ne; double *edge_coords = network->edge_coords; for (unsigned int j = 0; j < num_edges; j++) { if (edge_coords[2 * j + 1] < crack_len && fabs(edge_coords[2 * j] - network->L / 2) < crack_width) { instance->fuses[j] = true; instance->num_remaining_edges--; } else fuses[j] = false; } return true; }