From 2bb0740b68fdb62d45adc00204b3990ca42ade77 Mon Sep 17 00:00:00 2001 From: Jaron Kent-Dobias Date: Mon, 22 Aug 2016 10:11:14 -0400 Subject: started repo again without all the data files gunking the works --- src/update_boundary.c | 86 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 86 insertions(+) create mode 100644 src/update_boundary.c (limited to 'src/update_boundary.c') diff --git a/src/update_boundary.c b/src/update_boundary.c new file mode 100644 index 0000000..3d2a086 --- /dev/null +++ b/src/update_boundary.c @@ -0,0 +1,86 @@ + +#include "fracture.h" + +void update_boundary(finst *instance, const double *avg_field) { + int size = instance->network->num_edges; + int width = sqrt(size); + int num_verts = instance->network->num_verts; + double *boundary = (double *)instance->boundary_cond->x; + + boundary[num_verts] = 0; + boundary[num_verts + 1] = 0; + + if (instance->voltage_bound) { + for (int i = 0; i < width / 4; i++) { + int t_ind = (width + 1) * (width / 8) + width / 4 - 1 + i; + double t_val = avg_field[t_ind]; + boundary[2 * i] = t_val; + boundary[2 * i + 1] = t_val; + + int b_ind = num_verts - (width + 1) * (width / 8) - width / 4 - i; + double b_val = avg_field[b_ind]; + boundary[num_verts - 1 - 2 * i] = b_val; + boundary[num_verts - 2 - 2 * i] = b_val; + + int l_ind = (width + 1) * (width / 8) + width / 2 + width / 4 - 1 + + i * (width + 1); + double l_val = avg_field[l_ind]; + boundary[width / 2 + 2 * i * (width + 1)] = l_val; + boundary[width / 2 + (2 * i + 1) * (width + 1)] = l_val; + + int r_ind = (width + 1) * (width / 8) + width - 1 + i * (width + 1); + double r_val = avg_field[r_ind]; + boundary[width + 2 * i * (width + 1)] = r_val; + boundary[width + (2 * i + 1) * (width + 1)] = r_val; + } + } else { + const double *stress = avg_field; + + for (int i = 0; i < width / 4; i++) { + boundary[2 * i] = stress[(width / 4 - 1) * width + width / 4 + 2 * i] + + stress[(width / 4 - 1) * width + width / 4 + 1 + 2 * i]; + boundary[2 * i + 1] = + stress[(width / 4 - 1) * width + width / 4 + 2 * i] + + stress[(width / 4 - 1) * width + width / 4 + 1 + 2 * i]; + + boundary[num_verts - 2 * i - 1] = + -stress[size - 1 - ((width / 4 - 1) * width + width / 4 + 2 * i)] - + stress[size - 1 - ((width / 4 - 1) * width + width / 4 + 1 + 2 * i)]; + boundary[num_verts - 2 * i - 2] = + -stress[size - 1 - ((width / 4 - 1) * width + width / 4 + 2 * i)] - + stress[size - 1 - ((width / 4 - 1) * width + width / 4 + 1 + 2 * i)]; + + boundary[(width + 1) / 2 + 2 * i * (width + 1)] = + stress[((width / 4) + 2 * i) * width + width / 4 - 1] - + stress[((width / 4) + 2 * i + 1) * width + width / 4 - 1]; + boundary[(width + 1) / 2 + (2 * i + 1) * (width + 1)] = + stress[((width / 4) + 2 * i) * width + width / 4 - 1] - + stress[((width / 4) + 2 * i + 1) * width + width / 4 - 1]; + + boundary[(width + 1) / 2 + 2 * i * (width + 1) + width / 2] = + stress[((width / 4) + 2 * i) * width + width / 4 + width / 2] - + stress[((width / 4) + 2 * i + 1) * width + width / 4 + width / 2]; + boundary[(width + 1) / 2 + (2 * i + 1) * (width + 1) + width / 2] = + stress[((width / 4) + 2 * i) * width + width / 4 + width / 2] - + stress[((width / 4) + 2 * i + 1) * width + width / 4 + width / 2]; + } + + double s_total = 0; + double total = 0; + for (int i = 0; i < num_verts; i++) { + if (boundary[i] != boundary[i]) { + boundary[i] = 0; + } else { + total += fabs(boundary[i]); + s_total += boundary[i]; + } + } + + // assert(fabs(s_total) < inf); + boundary[num_verts] -= s_total; + + for (int i = 0; i < num_verts + 1; i++) { + boundary[i] /= total; + } + } +} -- cgit v1.2.3-70-g09d2