From 7ff906b9cd27a44472b40e78e5d595ea41df1482 Mon Sep 17 00:00:00 2001 From: pants Date: Wed, 31 Aug 2016 14:04:55 -0400 Subject: can generate voronoi networks with regular boundaries --- makefile_hal | 2 +- src/fortune/main.c | 6 ++-- src/fracture_network.c | 8 +++-- src/homo_voronoi_fracture.c | 33 ++++++++++-------- src/ini_network.c | 83 ++++++++++++++++++++++++++------------------- src/randfuncs.c | 22 ++++++------ 6 files changed, 88 insertions(+), 66 deletions(-) diff --git a/makefile_hal b/makefile_hal index 25ebf48..0b2f6ff 100644 --- a/makefile_hal +++ b/makefile_hal @@ -4,7 +4,7 @@ CFLAGS = -g -Os -O3 -Wall -fno-strict-aliasing -Wstrict-overflow -Wno-missing-fi LDFLAGS = -lcblas -llapack -ldl -lpthread -lcholmod -lamd -lcolamd -lsuitesparseconfig -lcamd -lccolamd -lm -lmetis -lgsl -lprofiler -ltcmalloc OBJ = break_data voronoi_bound_ini bin_values correlations beta_scales randfuncs instance get_dual_clusters coursegrain break_edge graph_components gen_laplacian geometry fracture_network get_current cracking_ini update_factor update_boundary get_file update_beta gen_voltcurmat ini_network free_network fortune/edgelist fortune/geometry fortune/heap fortune/main fortune/output fortune/voronoi fortune/memory get_conductivity -BIN = fracture current_scaling course_grain_square corr_test homo_voronoi_fracture homo_square_fracture compare_voronoi_fracture +BIN = current_scaling course_grain_square corr_test homo_voronoi_fracture homo_square_fracture compare_voronoi_fracture all: opt ${OBJ:%=obj/%.o} ${BIN:%=obj/%.o} ${BIN:%=bin/%} diff --git a/src/fortune/main.c b/src/fortune/main.c index 543f6e4..a1a4372 100644 --- a/src/fortune/main.c +++ b/src/fortune/main.c @@ -95,8 +95,8 @@ intptr_t *run_voronoi(unsigned int num, double *lattice, bool periodic, double x double *eff_lattice = lattice; if (periodic) { - unsigned int eff_num = 9 * num; - double *eff_lattice = (double *)malloc(2 * eff_num * sizeof(double)); + eff_num = 9 * num; + eff_lattice = (double *)malloc(2 * eff_num * sizeof(double)); for (unsigned int i = 0; i < num; i++) { // original sites - our baby boys @@ -139,7 +139,7 @@ intptr_t *run_voronoi(unsigned int num, double *lattice, bool periodic, double x nsites = eff_num; readsites(eff_lattice); - free(eff_lattice); + if (periodic) free(eff_lattice); next = nextone; siteidx = 0; diff --git a/src/fracture_network.c b/src/fracture_network.c index 148b08c..fe520fd 100644 --- a/src/fracture_network.c +++ b/src/fracture_network.c @@ -39,9 +39,11 @@ break_data *fracture_network(finst *instance, double *fuse_thres, } int last_broke = inc_break_fuses(instance, fuse_thres, field, cutoff); - if (last_broke > num_edges || last_broke < -1) { - printf("%g \n", conductivity); - getchar(); + if (last_broke > num_edges || last_broke < -1 || conductivity < 1e-8) { + printf("whoops %u\n\n", breaking_data->num_broken); + free(voltages); + free(field); + break; } update_break_data(breaking_data, last_broke, fabs(conductivity * fuse_thres[last_broke] / field[last_broke]), conductivity); diff --git a/src/homo_voronoi_fracture.c b/src/homo_voronoi_fracture.c index 65a25f8..26ceed4 100644 --- a/src/homo_voronoi_fracture.c +++ b/src/homo_voronoi_fracture.c @@ -7,7 +7,7 @@ int main(int argc, char *argv[]) { // defining variables to be (potentially) set by command line flags uint8_t filename_len; - uint64_t N; + uint32_t N; uint_t L; double beta, inf, cutoff; bool include_breaking, save_cluster_dist, use_voltage_boundaries, use_dual, save_network, @@ -63,6 +63,11 @@ int main(int argc, char *argv[]) { use_voltage_boundaries = true; boundc2 = 't'; break; + case 3: + boundary = EMBEDDED_BOUND; + boundc2 = 'e'; + use_dual = true; + break; default: printf("boundary specifier must be 0 (FREE_BOUND), 1 (CYLINDER_BOUND), or 2 (TORUS_BOUND).\n"); exit(EXIT_FAILURE); @@ -124,15 +129,15 @@ int main(int argc, char *argv[]) { } // define arrays for saving cluster and avalanche distributions - uint64_t *cluster_size_dist; - uint64_t *avalanche_size_dist; + uint32_t *cluster_size_dist; + uint32_t *avalanche_size_dist; char *c_filename; char *a_filename; if (save_cluster_dist) { cluster_size_dist = - (uint64_t *)malloc(c_dist_size * sizeof(uint64_t)); + (uint32_t *)malloc(c_dist_size * sizeof(uint32_t)); avalanche_size_dist = - (uint64_t *)malloc(a_dist_size * sizeof(uint64_t)); + (uint32_t *)malloc(a_dist_size * sizeof(uint32_t)); c_filename = (char *)malloc(filename_len * sizeof(char)); a_filename = (char *)malloc(filename_len * sizeof(char)); @@ -143,11 +148,11 @@ int main(int argc, char *argv[]) { FILE *avalanche_out = fopen(a_filename, "rb"); if (cluster_out != NULL) { - fread(cluster_size_dist, sizeof(uint64_t), c_dist_size, cluster_out); + fread(cluster_size_dist, sizeof(uint32_t), c_dist_size, cluster_out); fclose(cluster_out); } if (avalanche_out != NULL) { - fread(avalanche_size_dist, sizeof(uint64_t), a_dist_size, avalanche_out); + fread(avalanche_size_dist, sizeof(uint32_t), a_dist_size, avalanche_out); fclose(avalanche_out); } } @@ -163,11 +168,11 @@ int main(int argc, char *argv[]) { } // define arrays for saving damage distributions - uint64_t *damage; + uint32_t *damage; char *d_filename; if (save_damage) { damage = - (uint64_t *)malloc(a_dist_size * sizeof(uint64_t)); + (uint32_t *)malloc(a_dist_size * sizeof(uint32_t)); d_filename = (char *)malloc(filename_len * sizeof(char)); snprintf(d_filename, filename_len, "damg_v_%c_%c_%d_%g.dat", boundc, boundc2, L, beta); @@ -175,7 +180,7 @@ int main(int argc, char *argv[]) { FILE *damage_out = fopen(d_filename, "rb"); if (damage_out != NULL) { - fread(damage, sizeof(uint64_t), a_dist_size, damage_out); + fread(damage, sizeof(uint32_t), a_dist_size, damage_out); fclose(damage_out); } } @@ -203,7 +208,7 @@ int main(int argc, char *argv[]) { printf("\n"); - for (uint64_t i = 0; i < N; i++) { + for (uint32_t i = 0; i < N; i++) { printf("\033[F\033[JFRACTURE: %0*d / %d\n", (uint8_t)log10(N) + 1, i + 1, N); fnet *network = ini_voronoi_network(L, boundary, use_dual, genfunc_hyperuniform, &c); @@ -340,8 +345,8 @@ int main(int argc, char *argv[]) { FILE *cluster_out = fopen(c_filename, "wb"); FILE *avalanche_out = fopen(a_filename, "wb"); - fwrite(cluster_size_dist, sizeof(uint64_t), c_dist_size, cluster_out); - fwrite(avalanche_size_dist, sizeof(uint64_t), a_dist_size, avalanche_out); + fwrite(cluster_size_dist, sizeof(uint32_t), c_dist_size, cluster_out); + fwrite(avalanche_size_dist, sizeof(uint32_t), a_dist_size, avalanche_out); fclose(cluster_out); fclose(avalanche_out); @@ -374,7 +379,7 @@ int main(int argc, char *argv[]) { if (save_damage) { FILE *hdam_file = fopen(d_filename, "wb"); - fwrite(damage, sizeof(uint64_t), a_dist_size, hdam_file); + fwrite(damage, sizeof(uint32_t), a_dist_size, hdam_file); fclose(hdam_file); free(d_filename); free(damage); diff --git a/src/ini_network.c b/src/ini_network.c index edf1539..263b1b2 100644 --- a/src/ini_network.c +++ b/src/ini_network.c @@ -23,12 +23,14 @@ double *get_edge_coords(unsigned int num_edges, double *vert_coords, unsigned int *edges_to_verts) { double *output = (double *)malloc(2 * num_edges * sizeof(double)); -#pragma omp parallel for + #pragma omp parallel for for (unsigned int i = 0; i < num_edges; i++) { unsigned int v1, v2; double v1x, v1y, v2x, v2y, dx, dy; v1 = edges_to_verts[2 * i]; v2 = edges_to_verts[2 * i + 1]; + output[2 * i] = 0; + output[2 * i + 1] = 0; v1x = vert_coords[2 * v1]; v1y = vert_coords[2 * v1 + 1]; v2x = vert_coords[2 * v2]; @@ -270,7 +272,8 @@ unsigned int *get_voro_dual_edges(unsigned int num_edges, unsigned int *triangles) { unsigned int *dual_edges = (unsigned int *)malloc(2 * num_edges * sizeof(unsigned int)); -#pragma omp parallel for + unsigned int place = 0; + #pragma omp parallel for for (unsigned int i = 0; i < num_edges; i++) { unsigned int v1, v2; v1 = edges[2 * i]; @@ -285,8 +288,9 @@ unsigned int *get_voro_dual_edges(unsigned int num_edges, t21 = triangles[3 * v2 + k]; t22 = triangles[3 * v2 + ((k + 1) % 3)]; if ((t11 == t21 && t12 == t22) || (t11 == t22 && t12 == t21)) { - dual_edges[2 * i] = t11 < t12 ? t11 : t12; - dual_edges[2 * i + 1] = t11 < t12 ? t12 : t11; + dual_edges[2 * place] = t11 < t12 ? t11 : t12; + dual_edges[2 * place + 1] = t11 < t12 ? t12 : t11; + place++; found_match = true; break; } @@ -361,7 +365,7 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, // prune the edges of the lattice and assign boundary vertices based on the // desired boundary conditions unsigned int num_bounds; - unsigned num_verts; + unsigned int num_verts; double *vert_coords; unsigned int *bound_inds; unsigned int *bound_verts; @@ -379,9 +383,9 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, edges = (unsigned int *)malloc(2 * tmp_num_edges * sizeof(unsigned int)); dual_edges = (unsigned int *)malloc(2 * tmp_num_edges * sizeof(unsigned int)); unsigned int num_t, num_b, num_l, num_r; - bool *bound_t, *bound_b, *bound_l, *bound_r; + bool *bound_top, *bound_b, *bound_l, *bound_r; num_t = 0; num_b = 0; num_l = 0; num_r = 0; - bound_t = (bool *)calloc(num_verts, sizeof(bool)); + bound_top = (bool *)calloc(num_verts, sizeof(bool)); bound_b = (bool *)calloc(num_verts, sizeof(bool)); bound_l = (bool *)calloc(num_verts, sizeof(bool)); bound_r = (bool *)calloc(num_verts, sizeof(bool)); @@ -394,15 +398,15 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, dx = v1x - v2x; dy = v1y - v2y; if (fabs(dy) > 0.5) { if (dy > 0) { - if (!bound_t[v1] && !bound_l[v1] && !bound_r[v1]) { - bound_t[v1] = true; num_t++; + if (!bound_top[v1] && !bound_l[v1] && !bound_r[v1]) { + bound_top[v1] = true; num_t++; } if (!bound_b[v2] && !bound_l[v2] && !bound_r[v2]) { bound_b[v2] = true; num_b++; } } else { - if (!bound_t[v2] && !bound_l[v2] && !bound_r[v2]) { - bound_t[v2] = true; num_t++; + if (!bound_top[v2] && !bound_l[v2] && !bound_r[v2]) { + bound_top[v2] = true; num_t++; } if (!bound_b[v1] && !bound_l[v1] && !bound_r[v1]) { bound_b[v1] = true; num_b++; @@ -410,17 +414,17 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, } } else if (fabs(dx) > 0.5) { if (dx > 0) { - if (!bound_r[v1] && !bound_t[v1] && !bound_b[v1]) { + if (!bound_r[v1] && !bound_top[v1] && !bound_b[v1]) { bound_r[v1] = true; num_r++; } - if (!bound_l[v2] && !bound_t[v2] && !bound_b[v2]) { + if (!bound_l[v2] && !bound_top[v2] && !bound_b[v2]) { bound_l[v2] = true; num_l++; } } else { - if (!bound_r[v2] && !bound_t[v2] && !bound_b[v2]) { + if (!bound_r[v2] && !bound_top[v2] && !bound_b[v2]) { bound_r[v2] = true; num_r++; } - if (!bound_l[v1] && !bound_t[v1] && !bound_b[v1]) { + if (!bound_l[v1] && !bound_top[v1] && !bound_b[v1]) { bound_l[v1] = true; num_l++; } } @@ -441,7 +445,7 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, unsigned int pos_t, pos_b, pos_l, pos_r; pos_t = 0; pos_b = 0; pos_l = 0; pos_r = 0; for (unsigned int i = 0; i < num_verts; i++) { - if (bound_t[i]) { + if (bound_top[i]) { bound_verts[pos_t] = i; pos_t++; } else if (bound_b[i]) { bound_verts[num_t + pos_b] = i; pos_b++; @@ -451,9 +455,10 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, bound_verts[num_t + num_b + num_l + pos_r] = i; pos_r++; } } - free(bound_l); free(bound_r); free(bound_t); free(bound_b); + free(bound_l); free(bound_r); free(bound_top); free(bound_b); free(tmp_edges); free(tmp_dual_edges); + num_bounds = 2; network->edges_to_verts_break = edges; network->edges_to_verts = edges; network->num_verts_break = num_verts; @@ -471,9 +476,9 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, edges = (unsigned int *)malloc(2 * tmp_num_edges * sizeof(unsigned int)); dual_edges = (unsigned int *)malloc(2 * tmp_num_edges * sizeof(unsigned int)); unsigned int num_t, num_b; - bool *bound_t, *bound_b; + bool *bound_top, *bound_b; num_t = 0; num_b = 0; - bound_t = (bool *)calloc(num_verts, sizeof(bool)); + bound_top = (bool *)calloc(num_verts, sizeof(bool)); bound_b = (bool *)calloc(num_verts, sizeof(bool)); for (unsigned int i = 0; i < tmp_num_edges; i++) { unsigned int v1, v2; @@ -483,15 +488,15 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, dy = v1y - v2y; if (fabs(dy) > 0.5) { if (dy > 0) { - if (!bound_t[v1]) { - bound_t[v1] = true; num_t++; + if (!bound_top[v1]) { + bound_top[v1] = true; num_t++; } if (!bound_b[v2]) { bound_b[v2] = true; num_b++; } } else { - if (!bound_t[v2]) { - bound_t[v2] = true; num_t++; + if (!bound_top[v2]) { + bound_top[v2] = true; num_t++; } if (!bound_b[v1]) { bound_b[v1] = true; num_b++; @@ -512,13 +517,13 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, unsigned int pos_t, pos_b; pos_t = 0; pos_b = 0; for (unsigned int i = 0; i < num_verts; i++) { - if (bound_t[i]) { + if (bound_top[i]) { bound_verts[pos_t] = i; pos_t++; } else if (bound_b[i]) { bound_verts[num_t + pos_b] = i; pos_b++; } } - free(bound_t); free(bound_b); + free(bound_top); free(bound_b); free(tmp_edges); free(tmp_dual_edges); network->edges_to_verts_break = edges; @@ -539,7 +544,7 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, edges[2*i+1] = tmp_edges[2*i+1]; } dual_edges = tmp_dual_edges; - bool *bound_t = (bool *)calloc(tmp_num_verts, sizeof(bool)); + bool *bound_top = (bool *)calloc(tmp_num_verts, sizeof(bool)); int *edge_change = (int *)calloc(num_edges, sizeof(int)); unsigned int num_t = 0; for (unsigned int i = 0; i < num_edges; i++) { @@ -552,13 +557,13 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, if (fabs(dy) > 0.5) { if (dy > 0) { edge_change[i] = 1; - if (!bound_t[v1]) { - bound_t[v1] = true; num_t++; + if (!bound_top[v1]) { + bound_top[v1] = true; num_t++; } } else { edge_change[i] = 2; - if (!bound_t[v2]) { - bound_t[v2] = true; num_t++; + if (!bound_top[v2]) { + bound_top[v2] = true; num_t++; } } } @@ -571,7 +576,7 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, for (unsigned int i = 0; i < tmp_num_verts; i++) { vert_coords[2*i] = tmp_vert_coords[2*i]; vert_coords[2*i+1] = tmp_vert_coords[2*i+1]; - if (bound_t[i]) { + if (bound_top[i]) { bound_verts[pos_t] = i; vert_coords[2*(tmp_num_verts + pos_t)] = tmp_vert_coords[2*i]; vert_coords[2*(tmp_num_verts + pos_t)+1] = tmp_vert_coords[2*i+1]; @@ -589,7 +594,7 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, } } free(tmp_vert_coords); - free(bound_t); + free(bound_top); free(edge_change); network->num_verts_break = num_verts; network->num_verts = tmp_num_verts; @@ -599,17 +604,25 @@ fnet *ini_voronoi_network(unsigned int L, bound_t boundary, bool use_dual, break; } case EMBEDDED_BOUND: { - num_bounds = 4; + num_bounds = 2; bound_inds = (unsigned int *)malloc(5 * sizeof(unsigned int)); bound_verts = (unsigned int *)malloc(2 * L * sizeof(unsigned int)); for (unsigned int i = 0; i < 5; i++) bound_inds[i] = i * L / 2; for (unsigned int i = 0; i < 2 * L; i++) bound_verts[i] = i; - num_edges = tmp_num_edges; + unsigned int num_away = 0; + for (unsigned int i = 0; i < tmp_num_edges; i++) { + if (tmp_dual_edges[2*i] > num || tmp_dual_edges[2*i+1] > num) num_away++; + } + num_edges = (int)tmp_num_edges - (int)num_away; num_verts = tmp_num_verts; edges = tmp_edges; dual_edges = tmp_dual_edges; - num_edges = tmp_num_edges; vert_coords = tmp_vert_coords; + network->num_verts_break = num_verts; + network->num_verts = num_verts; + network->edges_to_verts_break = edges; + network->edges_to_verts = edges; + network->break_dim = num_verts + num_bounds; } } diff --git a/src/randfuncs.c b/src/randfuncs.c index 35b9f56..830e078 100644 --- a/src/randfuncs.c +++ b/src/randfuncs.c @@ -22,23 +22,25 @@ double *genfunc_hyperuniform(unsigned int L, bound_t boundary, gsl_rng *r, unsig double *lattice = (double *)malloc(2 * (*num) * sizeof(double)); double rho = *num; unsigned int to_gen = *num; + unsigned int start = 0; if (boundary == EMBEDDED_BOUND) { for (unsigned int i = 0; i < L / 2; i++) { - lattice[2 * i] = 0; - lattice[2 * i + 1] = (2. * i + 1.) / L; - - lattice[L / 2 + 2 * i] = 1; - lattice[L / 2 + 2 * i + 1] = (2. * i + 1.) / L; + lattice[2 * i + 1] = 0; + lattice[2 * i] = (2. * i + 1.) / L; + lattice[L + 2 * i + 1] = 1; lattice[L + 2 * i] = (2. * i + 1.) / L; - lattice[L + 2 * i + 1] = 0; - lattice[3 * L / 2 + 2 * i] = (2. * i + 1.) / L; - lattice[3 * L / 2 + 2 * i + 1] = 1; + lattice[2 * L + 2 * i + 1] = (2. * i + 1.) / L; + lattice[2 * L + 2 * i] = 0; + + lattice[3 * L + 2 * i + 1] = (2. * i + 1.) / L; + lattice[3 * L + 2 * i] = 1; } to_gen -= 2 * L; + start = 2 * L; } for (unsigned int i = 0; i < to_gen; i++) { @@ -67,8 +69,8 @@ double *genfunc_hyperuniform(unsigned int L, bound_t boundary, gsl_rng *r, unsig } } } - lattice[2 * i] = x; - lattice[2 * i + 1] = y; + lattice[2*start + 2 * i] = x; + lattice[2*start + 2 * i + 1] = y; } return lattice; -- cgit v1.2.3-70-g09d2