diff options
| author | Jaron Kent-Dobias <jaron@kent-dobias.com> | 2017-06-07 16:15:03 -0400 | 
|---|---|---|
| committer | Jaron Kent-Dobias <jaron@kent-dobias.com> | 2017-06-07 16:15:03 -0400 | 
| commit | fa590b7542c11b417ee87eebdd162b7aa08b57d4 (patch) | |
| tree | d4b884171ccd1072a78cfe0229b196523760b2c8 /src | |
| parent | 23682e445ba64452c453346f516f9101761438ba (diff) | |
| download | fuse_networks-fa590b7542c11b417ee87eebdd162b7aa08b57d4.tar.gz fuse_networks-fa590b7542c11b417ee87eebdd162b7aa08b57d4.tar.bz2 fuse_networks-fa590b7542c11b417ee87eebdd162b7aa08b57d4.zip  | |
ran clang-format on all files
Diffstat (limited to 'src')
| -rw-r--r-- | src/anal_process.c | 215 | ||||
| -rw-r--r-- | src/big_anal_process.c | 244 | ||||
| -rw-r--r-- | src/cen_anal_process.c | 245 | ||||
| -rw-r--r-- | src/corr_test.c | 30 | ||||
| -rw-r--r-- | src/fracture.c | 891 | ||||
| -rw-r--r-- | src/fracture.h | 73 | ||||
| -rw-r--r-- | src/long_anal_process.c | 244 | 
7 files changed, 978 insertions, 964 deletions
diff --git a/src/anal_process.c b/src/anal_process.c index 0ebec84..de27571 100644 --- a/src/anal_process.c +++ b/src/anal_process.c @@ -1,134 +1,135 @@  #include "fracture.h" +#include <gsl/gsl_blas.h> +#include <gsl/gsl_matrix.h>  #include <gsl/gsl_sf_erf.h>  #include <gsl/gsl_sf_laguerre.h> -#include <gsl/gsl_matrix.h>  #include <gsl/gsl_vector.h> -#include <gsl/gsl_blas.h>  void mom(uint_t len, uint_t n, uint32_t *data, double result[2]) { -	uint_t total = 0; -	double mom = 0; -	double mom_err = 0; +  uint_t total = 0; +  double mom = 0; +  double mom_err = 0; -	for (uint_t i = 0; i < len; i++) { -		uint32_t datai = data[i]; -		double in = pow(i, n); +  for (uint_t i = 0; i < len; i++) { +    uint32_t datai = data[i]; +    double in = pow(i, n); -		total += datai; -		mom += in * datai; -		mom_err += gsl_pow_2(in) * datai;  -	} +    total += datai; +    mom += in * datai; +    mom_err += gsl_pow_2(in) * datai; +  } -	double momf = mom / total; -	double momf_err = momf * sqrt(mom_err / gsl_pow_2(mom) + 1 / total); +  double momf = mom / total; +  double momf_err = momf * sqrt(mom_err / gsl_pow_2(mom) + 1 / total); -	result[0] = momf; -	result[1] = momf_err; +  result[0] = momf; +  result[1] = momf_err;  }  int main(int argc, char *argv[]) { -	uint_t nc = argc - 1; -	uint_t *Ls_c = (uint_t *)malloc(nc * sizeof(uint_t)); -	double *betas_c = (double *)malloc(nc * sizeof(double)); -	double *vals_c1 = (double *)malloc(nc * sizeof(double)); -	double *errors_c1 = (double *)malloc(nc * sizeof(double)); -	double *vals_c2 = (double *)malloc(nc * sizeof(double)); -	double *errors_c2 = (double *)malloc(nc * sizeof(double)); -	double *vals_c3 = (double *)malloc(nc * sizeof(double)); -	double *errors_c3 = (double *)malloc(nc * sizeof(double)); -	double *vals_c4 = (double *)malloc(nc * sizeof(double)); -	double *errors_c4 = (double *)malloc(nc * sizeof(double)); - -	char *out_filename = (char *)malloc(100 * sizeof(char)); -	uint_t out_filename_len = 0; +  uint_t nc = argc - 1; +  uint_t *Ls_c = (uint_t *)malloc(nc * sizeof(uint_t)); +  double *betas_c = (double *)malloc(nc * sizeof(double)); +  double *vals_c1 = (double *)malloc(nc * sizeof(double)); +  double *errors_c1 = (double *)malloc(nc * sizeof(double)); +  double *vals_c2 = (double *)malloc(nc * sizeof(double)); +  double *errors_c2 = (double *)malloc(nc * sizeof(double)); +  double *vals_c3 = (double *)malloc(nc * sizeof(double)); +  double *errors_c3 = (double *)malloc(nc * sizeof(double)); +  double *vals_c4 = (double *)malloc(nc * sizeof(double)); +  double *errors_c4 = (double *)malloc(nc * sizeof(double)); -	for (uint_t i = 0; i < nc; i++) { -		char *fn = argv[1 + i]; -		char *buff = (char *)malloc(20 * sizeof(char)); -		uint_t pos = 0; uint_t c = 0; -		while (c < 5) { -			if (fn[pos] == '_') { -				c++; -			} -			if (i == 0) { -				out_filename[pos] = fn[pos]; -				out_filename_len++; -			} -			pos++; -		} -		c = 0; -		while (fn[pos] != '_') { -			buff[c] = fn[pos]; -			pos++; -			c++; -		} -		buff[c] = '\0'; -		Ls_c[i] = atoi(buff); -		c = 0; -		pos++; -		while (fn[pos] != '_') { -			buff[c] = fn[pos]; -			pos++; -			c++; -		} -		buff[c] = '\0'; -		betas_c[i] = atof(buff); -		free(buff); +  char *out_filename = (char *)malloc(100 * sizeof(char)); +  uint_t out_filename_len = 0; -		uint_t dist_len = 4 * pow(Ls_c[i], 2); -		uint32_t *dist = malloc(dist_len * sizeof(uint32_t)); -		FILE *dist_file = fopen(fn, "rb"); -		fread(dist, sizeof(uint32_t), dist_len, dist_file); -		fclose(dist_file); -		{ -			double mom1[2]; -			mom(dist_len, 1, dist, mom1); -			vals_c1[i] = mom1[0]; -			errors_c1[i] = mom1[1]; +  for (uint_t i = 0; i < nc; i++) { +    char *fn = argv[1 + i]; +    char *buff = (char *)malloc(20 * sizeof(char)); +    uint_t pos = 0; +    uint_t c = 0; +    while (c < 5) { +      if (fn[pos] == '_') { +        c++; +      } +      if (i == 0) { +        out_filename[pos] = fn[pos]; +        out_filename_len++; +      } +      pos++; +    } +    c = 0; +    while (fn[pos] != '_') { +      buff[c] = fn[pos]; +      pos++; +      c++; +    } +    buff[c] = '\0'; +    Ls_c[i] = atoi(buff); +    c = 0; +    pos++; +    while (fn[pos] != '_') { +      buff[c] = fn[pos]; +      pos++; +      c++; +    } +    buff[c] = '\0'; +    betas_c[i] = atof(buff); +    free(buff); -			double mom2[2]; -			mom(dist_len, 2, dist, mom2); -			vals_c2[i] = mom2[0]; -			errors_c2[i] = mom2[1]; +    uint_t dist_len = 4 * pow(Ls_c[i], 2); +    uint32_t *dist = malloc(dist_len * sizeof(uint32_t)); +    FILE *dist_file = fopen(fn, "rb"); +    fread(dist, sizeof(uint32_t), dist_len, dist_file); +    fclose(dist_file); +    { +      double mom1[2]; +      mom(dist_len, 1, dist, mom1); +      vals_c1[i] = mom1[0]; +      errors_c1[i] = mom1[1]; -			double mom3[2]; -			mom(dist_len, 3, dist, mom3); -			vals_c3[i] = mom3[0]; -			errors_c3[i] = mom3[1]; +      double mom2[2]; +      mom(dist_len, 2, dist, mom2); +      vals_c2[i] = mom2[0]; +      errors_c2[i] = mom2[1]; -			double mom4[2]; -			mom(dist_len, 4, dist, mom4); -			vals_c4[i] = mom4[0]; -			errors_c4[i] = mom4[1]; -		} -		free(dist); -	} +      double mom3[2]; +      mom(dist_len, 3, dist, mom3); +      vals_c3[i] = mom3[0]; +      errors_c3[i] = mom3[1]; -	out_filename[out_filename_len-1] = '.'; -	out_filename[out_filename_len] = 't'; -	out_filename[out_filename_len+1] = 'x'; -	out_filename[out_filename_len+2] = 't'; -	out_filename[out_filename_len+3] = '\0'; +      double mom4[2]; +      mom(dist_len, 4, dist, mom4); +      vals_c4[i] = mom4[0]; +      errors_c4[i] = mom4[1]; +    } +    free(dist); +  } -	FILE *out_file = fopen(out_filename, "w"); +  out_filename[out_filename_len - 1] = '.'; +  out_filename[out_filename_len] = 't'; +  out_filename[out_filename_len + 1] = 'x'; +  out_filename[out_filename_len + 2] = 't'; +  out_filename[out_filename_len + 3] = '\0'; -	for (uint_t i = 0; i < nc; i++) { -		fprintf(out_file, "%u %g %g %g %g %g %g %g %g %g\n", Ls_c[i], betas_c[i], vals_c1[i], errors_c1[i], vals_c2[i], errors_c2[i], vals_c3[i], errors_c3[i], vals_c4[i], errors_c4[i]); -	} +  FILE *out_file = fopen(out_filename, "w"); -	fclose(out_file); +  for (uint_t i = 0; i < nc; i++) { +    fprintf(out_file, "%u %g %g %g %g %g %g %g %g %g\n", Ls_c[i], betas_c[i], +            vals_c1[i], errors_c1[i], vals_c2[i], errors_c2[i], vals_c3[i], +            errors_c3[i], vals_c4[i], errors_c4[i]); +  } +  fclose(out_file); -	free(Ls_c); -	free(betas_c); -	free(vals_c1); -	free(errors_c1); -	free(vals_c2); -	free(errors_c2); -	free(vals_c3); -	free(errors_c3); +  free(Ls_c); +  free(betas_c); +  free(vals_c1); +  free(errors_c1); +  free(vals_c2); +  free(errors_c2); +  free(vals_c3); +  free(errors_c3); -	return 0; +  return 0;  } - diff --git a/src/big_anal_process.c b/src/big_anal_process.c index 0f7106f..8c1d8ba 100644 --- a/src/big_anal_process.c +++ b/src/big_anal_process.c @@ -1,156 +1,158 @@  #include "fracture.h" +#include <gsl/gsl_blas.h> +#include <gsl/gsl_matrix.h>  #include <gsl/gsl_sf_erf.h>  #include <gsl/gsl_sf_laguerre.h> -#include <gsl/gsl_matrix.h>  #include <gsl/gsl_vector.h> -#include <gsl/gsl_blas.h>  #include <sys/stat.h>  void get_mean(uint_t len, double *data, double result[2]) { -	double total = 0; +  double total = 0; -	for (uint_t i = 0; i < len; i++) { -		total += data[i]; -	} +  for (uint_t i = 0; i < len; i++) { +    total += data[i]; +  } -	double mean = total / len; -	double total_sq_diff = 0; +  double mean = total / len; +  double total_sq_diff = 0; -	for (uint_t i = 0; i < len; i++) { -		total_sq_diff += pow(mean - data[i], 2); -	} +  for (uint_t i = 0; i < len; i++) { +    total_sq_diff += pow(mean - data[i], 2); +  } -	double mean_err = sqrt(total_sq_diff) / len; +  double mean_err = sqrt(total_sq_diff) / len; -	result[0] = mean; -	result[1] = mean_err; +  result[0] = mean; +  result[1] = mean_err;  } -void get_mom(uint_t len, uint_t n, double *data, double mean[2], double result[2]) { -	double total_n_diff = 0; -	double total_n2_diff = 0; +void get_mom(uint_t len, uint_t n, double *data, double mean[2], +             double result[2]) { +  double total_n_diff = 0; +  double total_n2_diff = 0; -	for (uint_t i = 0; i < len; i++) { -		total_n_diff += pow(fabs(mean[0] - data[i]), n); -		total_n2_diff += pow(fabs(mean[0] - data[i]), 2 * n); -	} +  for (uint_t i = 0; i < len; i++) { +    total_n_diff += pow(fabs(mean[0] - data[i]), n); +    total_n2_diff += pow(fabs(mean[0] - data[i]), 2 * n); +  } -	double mom = total_n_diff / len; -	double mom_err = sqrt(total_n2_diff) / len; +  double mom = total_n_diff / len; +  double mom_err = sqrt(total_n2_diff) / len; -	result[0] = mom; -	result[1] = mom_err; +  result[0] = mom; +  result[1] = mom_err;  }  int main(int argc, char *argv[]) { -	uint_t nc = argc - 1; -	uint_t *Ls_c = (uint_t *)malloc(nc * sizeof(uint_t)); -	double *betas_c = (double *)malloc(nc * sizeof(double)); -	double *vals_c1 = (double *)malloc(nc * sizeof(double)); -	double *errors_c1 = (double *)malloc(nc * sizeof(double)); -	double *vals_c2 = (double *)malloc(nc * sizeof(double)); -	double *errors_c2 = (double *)malloc(nc * sizeof(double)); -	double *vals_c3 = (double *)malloc(nc * sizeof(double)); -	double *errors_c3 = (double *)malloc(nc * sizeof(double)); -	double *vals_c4 = (double *)malloc(nc * sizeof(double)); -	double *errors_c4 = (double *)malloc(nc * sizeof(double)); - -	char *out_filename = (char *)malloc(100 * sizeof(char)); -	uint_t out_filename_len = 0; +  uint_t nc = argc - 1; +  uint_t *Ls_c = (uint_t *)malloc(nc * sizeof(uint_t)); +  double *betas_c = (double *)malloc(nc * sizeof(double)); +  double *vals_c1 = (double *)malloc(nc * sizeof(double)); +  double *errors_c1 = (double *)malloc(nc * sizeof(double)); +  double *vals_c2 = (double *)malloc(nc * sizeof(double)); +  double *errors_c2 = (double *)malloc(nc * sizeof(double)); +  double *vals_c3 = (double *)malloc(nc * sizeof(double)); +  double *errors_c3 = (double *)malloc(nc * sizeof(double)); +  double *vals_c4 = (double *)malloc(nc * sizeof(double)); +  double *errors_c4 = (double *)malloc(nc * sizeof(double)); -	for (uint_t i = 0; i < nc; i++) { -		char *fn = argv[1 + i]; -		char *buff = (char *)malloc(20 * sizeof(char)); -		uint_t pos = 0; uint_t c = 0; -		while (c < 5) { -			if (fn[pos] == '_') { -				c++; -			} -			if (i == 0) { -				out_filename[pos] = fn[pos]; -				out_filename_len++; -			} -			pos++; -		} -		c = 0; -		while (fn[pos] != '_') { -			buff[c] = fn[pos]; -			pos++; -			c++; -		} -		buff[c] = '\0'; -		Ls_c[i] = atoi(buff); -		c = 0; -		pos++; -		while (fn[pos] != '_') { -			buff[c] = fn[pos]; -			pos++; -			c++; -		} -		buff[c] = '\0'; -		betas_c[i] = atof(buff); -		free(buff); +  char *out_filename = (char *)malloc(100 * sizeof(char)); +  uint_t out_filename_len = 0; -		struct stat info; -		stat(fn, &info); -		uint_t num_bytes = info.st_size; -		uint_t dist_len = (num_bytes * sizeof(char)) / sizeof(double); +  for (uint_t i = 0; i < nc; i++) { +    char *fn = argv[1 + i]; +    char *buff = (char *)malloc(20 * sizeof(char)); +    uint_t pos = 0; +    uint_t c = 0; +    while (c < 5) { +      if (fn[pos] == '_') { +        c++; +      } +      if (i == 0) { +        out_filename[pos] = fn[pos]; +        out_filename_len++; +      } +      pos++; +    } +    c = 0; +    while (fn[pos] != '_') { +      buff[c] = fn[pos]; +      pos++; +      c++; +    } +    buff[c] = '\0'; +    Ls_c[i] = atoi(buff); +    c = 0; +    pos++; +    while (fn[pos] != '_') { +      buff[c] = fn[pos]; +      pos++; +      c++; +    } +    buff[c] = '\0'; +    betas_c[i] = atof(buff); +    free(buff); -		double *dist = malloc(dist_len * sizeof(double)); -		FILE *dist_file = fopen(fn, "rb"); -		fread(dist, sizeof(double), dist_len, dist_file); -		fclose(dist_file); -		{ -			double mom1[2]; -			get_mean(dist_len, dist, mom1); -			vals_c1[i] = mom1[0]; -			errors_c1[i] = mom1[1]; +    struct stat info; +    stat(fn, &info); +    uint_t num_bytes = info.st_size; +    uint_t dist_len = (num_bytes * sizeof(char)) / sizeof(double); -			double mom2[2]; -			get_mom(dist_len, 2, dist, mom1, mom2); -			vals_c2[i] = mom2[0]; -			errors_c2[i] = mom2[1]; +    double *dist = malloc(dist_len * sizeof(double)); +    FILE *dist_file = fopen(fn, "rb"); +    fread(dist, sizeof(double), dist_len, dist_file); +    fclose(dist_file); +    { +      double mom1[2]; +      get_mean(dist_len, dist, mom1); +      vals_c1[i] = mom1[0]; +      errors_c1[i] = mom1[1]; -			double mom3[2]; -			get_mom(dist_len, 3, dist, mom1, mom3); -			vals_c3[i] = mom3[0]; -			errors_c3[i] = mom3[1]; +      double mom2[2]; +      get_mom(dist_len, 2, dist, mom1, mom2); +      vals_c2[i] = mom2[0]; +      errors_c2[i] = mom2[1]; -			double mom4[2]; -			get_mom(dist_len, 4, dist, mom1, mom4); -			vals_c4[i] = mom4[0]; -			errors_c4[i] = mom4[1]; -		} -		free(dist); -	} +      double mom3[2]; +      get_mom(dist_len, 3, dist, mom1, mom3); +      vals_c3[i] = mom3[0]; +      errors_c3[i] = mom3[1]; -	out_filename[out_filename_len-1] = '.'; -	out_filename[out_filename_len] = 't'; -	out_filename[out_filename_len+1] = 'x'; -	out_filename[out_filename_len+2] = 't'; -	out_filename[out_filename_len+3] = '\0'; +      double mom4[2]; +      get_mom(dist_len, 4, dist, mom1, mom4); +      vals_c4[i] = mom4[0]; +      errors_c4[i] = mom4[1]; +    } +    free(dist); +  } -	FILE *out_file = fopen(out_filename, "w"); +  out_filename[out_filename_len - 1] = '.'; +  out_filename[out_filename_len] = 't'; +  out_filename[out_filename_len + 1] = 'x'; +  out_filename[out_filename_len + 2] = 't'; +  out_filename[out_filename_len + 3] = '\0'; -	for (uint_t i = 0; i < nc; i++) { -		fprintf(out_file, "%u %g %g %g %g %g %g %g %g %g\n", Ls_c[i], betas_c[i], vals_c1[i], errors_c1[i], vals_c2[i], errors_c2[i], vals_c3[i], errors_c3[i], vals_c4[i], errors_c4[i]); -	} +  FILE *out_file = fopen(out_filename, "w"); -	fclose(out_file); +  for (uint_t i = 0; i < nc; i++) { +    fprintf(out_file, "%u %g %g %g %g %g %g %g %g %g\n", Ls_c[i], betas_c[i], +            vals_c1[i], errors_c1[i], vals_c2[i], errors_c2[i], vals_c3[i], +            errors_c3[i], vals_c4[i], errors_c4[i]); +  } +  fclose(out_file); -	free(Ls_c); -	free(betas_c); -	free(vals_c1); -	free(errors_c1); -	free(vals_c2); -	free(errors_c2); -	free(vals_c3); -	free(errors_c3); -	free(vals_c4); -	free(errors_c4); +  free(Ls_c); +  free(betas_c); +  free(vals_c1); +  free(errors_c1); +  free(vals_c2); +  free(errors_c2); +  free(vals_c3); +  free(errors_c3); +  free(vals_c4); +  free(errors_c4); -	return 0; +  return 0;  } - diff --git a/src/cen_anal_process.c b/src/cen_anal_process.c index 3bf388c..ee2b029 100644 --- a/src/cen_anal_process.c +++ b/src/cen_anal_process.c @@ -1,154 +1,157 @@  #include "fracture.h" +#include <gsl/gsl_blas.h> +#include <gsl/gsl_matrix.h>  #include <gsl/gsl_sf_erf.h>  #include <gsl/gsl_sf_laguerre.h> -#include <gsl/gsl_matrix.h>  #include <gsl/gsl_vector.h> -#include <gsl/gsl_blas.h>  void get_mean(uint_t len, uint32_t *data, double result[2]) { -	uint_t total = 0; -	double mean = 0; -	double mean_err = 0; +  uint_t total = 0; +  double mean = 0; +  double mean_err = 0; -	for (uint_t i = 0; i < len; i++) { -		uint32_t datai = data[i]; +  for (uint_t i = 0; i < len; i++) { +    uint32_t datai = data[i]; -		total += datai; -		mean += i * datai; -		mean_err += gsl_pow_2(i) * datai;  -	} +    total += datai; +    mean += i * datai; +    mean_err += gsl_pow_2(i) * datai; +  } -	double meanf = mean / total; -	double meanf_err = meanf * sqrt(mean_err / gsl_pow_2(mean) + 1 / total); +  double meanf = mean / total; +  double meanf_err = meanf * sqrt(mean_err / gsl_pow_2(mean) + 1 / total); -	result[0] = meanf; -	result[1] = meanf_err; +  result[0] = meanf; +  result[1] = meanf_err;  } -void get_mom(uint_t len, uint_t n, uint32_t *data, double mean[2], double result[2]) { -	uint_t total = 0; -	double mom = 0; -	double mom_err = 0; +void get_mom(uint_t len, uint_t n, uint32_t *data, double mean[2], +             double result[2]) { +  uint_t total = 0; +  double mom = 0; +  double mom_err = 0; -	for (uint_t i = 0; i < len; i++) { -		uint32_t datai = data[i]; -		double in = pow(fabs(((double)i) - mean[0]), n); +  for (uint_t i = 0; i < len; i++) { +    uint32_t datai = data[i]; +    double in = pow(fabs(((double)i) - mean[0]), n); -		total += datai; -		mom += in * datai; -		mom_err += gsl_pow_2(in) * datai; // + gsl_pow_2(n * mean[1] / (((double)i) - mean[0]))); -	} +    total += datai; +    mom += in * datai; +    mom_err += gsl_pow_2(in) * +               datai; // + gsl_pow_2(n * mean[1] / (((double)i) - mean[0]))); +  } -	double momf = mom / total; -	double momf_err = momf * sqrt(mom_err / gsl_pow_2(mom) + 1 / total); +  double momf = mom / total; +  double momf_err = momf * sqrt(mom_err / gsl_pow_2(mom) + 1 / total); -	result[0] = momf; -	result[1] = momf_err; +  result[0] = momf; +  result[1] = momf_err;  }  int main(int argc, char *argv[]) { -	uint_t nc = argc - 1; -	uint_t *Ls_c = (uint_t *)malloc(nc * sizeof(uint_t)); -	double *betas_c = (double *)malloc(nc * sizeof(double)); -	double *vals_c1 = (double *)malloc(nc * sizeof(double)); -	double *errors_c1 = (double *)malloc(nc * sizeof(double)); -	double *vals_c2 = (double *)malloc(nc * sizeof(double)); -	double *errors_c2 = (double *)malloc(nc * sizeof(double)); -	double *vals_c3 = (double *)malloc(nc * sizeof(double)); -	double *errors_c3 = (double *)malloc(nc * sizeof(double)); -	double *vals_c4 = (double *)malloc(nc * sizeof(double)); -	double *errors_c4 = (double *)malloc(nc * sizeof(double)); - -	char *out_filename = (char *)malloc(100 * sizeof(char)); -	uint_t out_filename_len = 0; +  uint_t nc = argc - 1; +  uint_t *Ls_c = (uint_t *)malloc(nc * sizeof(uint_t)); +  double *betas_c = (double *)malloc(nc * sizeof(double)); +  double *vals_c1 = (double *)malloc(nc * sizeof(double)); +  double *errors_c1 = (double *)malloc(nc * sizeof(double)); +  double *vals_c2 = (double *)malloc(nc * sizeof(double)); +  double *errors_c2 = (double *)malloc(nc * sizeof(double)); +  double *vals_c3 = (double *)malloc(nc * sizeof(double)); +  double *errors_c3 = (double *)malloc(nc * sizeof(double)); +  double *vals_c4 = (double *)malloc(nc * sizeof(double)); +  double *errors_c4 = (double *)malloc(nc * sizeof(double)); -	for (uint_t i = 0; i < nc; i++) { -		char *fn = argv[1 + i]; -		char *buff = (char *)malloc(20 * sizeof(char)); -		uint_t pos = 0; uint_t c = 0; -		while (c < 5) { -			if (fn[pos] == '_') { -				c++; -			} -			if (i == 0) { -				out_filename[pos] = fn[pos]; -				out_filename_len++; -			} -			pos++; -		} -		c = 0; -		while (fn[pos] != '_') { -			buff[c] = fn[pos]; -			pos++; -			c++; -		} -		buff[c] = '\0'; -		Ls_c[i] = atoi(buff); -		c = 0; -		pos++; -		while (fn[pos] != '_') { -			buff[c] = fn[pos]; -			pos++; -			c++; -		} -		buff[c] = '\0'; -		betas_c[i] = atof(buff); -		free(buff); +  char *out_filename = (char *)malloc(100 * sizeof(char)); +  uint_t out_filename_len = 0; -		uint_t dist_len = 4 * pow(Ls_c[i], 2); -		uint32_t *dist = malloc(dist_len * sizeof(uint32_t)); -		FILE *dist_file = fopen(fn, "rb"); -		fread(dist, sizeof(uint32_t), dist_len, dist_file); -		fclose(dist_file); -		{ -			double mom1[2]; -			get_mean(dist_len, dist, mom1); -			vals_c1[i] = mom1[0]; -			errors_c1[i] = mom1[1]; +  for (uint_t i = 0; i < nc; i++) { +    char *fn = argv[1 + i]; +    char *buff = (char *)malloc(20 * sizeof(char)); +    uint_t pos = 0; +    uint_t c = 0; +    while (c < 5) { +      if (fn[pos] == '_') { +        c++; +      } +      if (i == 0) { +        out_filename[pos] = fn[pos]; +        out_filename_len++; +      } +      pos++; +    } +    c = 0; +    while (fn[pos] != '_') { +      buff[c] = fn[pos]; +      pos++; +      c++; +    } +    buff[c] = '\0'; +    Ls_c[i] = atoi(buff); +    c = 0; +    pos++; +    while (fn[pos] != '_') { +      buff[c] = fn[pos]; +      pos++; +      c++; +    } +    buff[c] = '\0'; +    betas_c[i] = atof(buff); +    free(buff); -			double mom2[2]; -			get_mom(dist_len, 2, dist, mom1, mom2); -			vals_c2[i] = mom2[0]; -			errors_c2[i] = mom2[1]; +    uint_t dist_len = 4 * pow(Ls_c[i], 2); +    uint32_t *dist = malloc(dist_len * sizeof(uint32_t)); +    FILE *dist_file = fopen(fn, "rb"); +    fread(dist, sizeof(uint32_t), dist_len, dist_file); +    fclose(dist_file); +    { +      double mom1[2]; +      get_mean(dist_len, dist, mom1); +      vals_c1[i] = mom1[0]; +      errors_c1[i] = mom1[1]; -			double mom3[2]; -			get_mom(dist_len, 3, dist, mom1, mom3); -			vals_c3[i] = mom3[0]; -			errors_c3[i] = mom3[1]; +      double mom2[2]; +      get_mom(dist_len, 2, dist, mom1, mom2); +      vals_c2[i] = mom2[0]; +      errors_c2[i] = mom2[1]; -			double mom4[2]; -			get_mom(dist_len, 4, dist, mom1, mom4); -			vals_c4[i] = mom4[0]; -			errors_c4[i] = mom4[1]; -		} -		free(dist); -	} +      double mom3[2]; +      get_mom(dist_len, 3, dist, mom1, mom3); +      vals_c3[i] = mom3[0]; +      errors_c3[i] = mom3[1]; -	out_filename[out_filename_len-1] = '.'; -	out_filename[out_filename_len] = 't'; -	out_filename[out_filename_len+1] = 'x'; -	out_filename[out_filename_len+2] = 't'; -	out_filename[out_filename_len+3] = '\0'; +      double mom4[2]; +      get_mom(dist_len, 4, dist, mom1, mom4); +      vals_c4[i] = mom4[0]; +      errors_c4[i] = mom4[1]; +    } +    free(dist); +  } -	FILE *out_file = fopen(out_filename, "w"); +  out_filename[out_filename_len - 1] = '.'; +  out_filename[out_filename_len] = 't'; +  out_filename[out_filename_len + 1] = 'x'; +  out_filename[out_filename_len + 2] = 't'; +  out_filename[out_filename_len + 3] = '\0'; -	for (uint_t i = 0; i < nc; i++) { -		fprintf(out_file, "%u %g %g %g %g %g %g %g %g %g\n", Ls_c[i], betas_c[i], vals_c1[i], errors_c1[i], vals_c2[i], errors_c2[i], vals_c3[i], errors_c3[i], vals_c4[i], errors_c4[i]); -	} +  FILE *out_file = fopen(out_filename, "w"); -	fclose(out_file); +  for (uint_t i = 0; i < nc; i++) { +    fprintf(out_file, "%u %g %g %g %g %g %g %g %g %g\n", Ls_c[i], betas_c[i], +            vals_c1[i], errors_c1[i], vals_c2[i], errors_c2[i], vals_c3[i], +            errors_c3[i], vals_c4[i], errors_c4[i]); +  } +  fclose(out_file); -	free(Ls_c); -	free(betas_c); -	free(vals_c1); -	free(errors_c1); -	free(vals_c2); -	free(errors_c2); -	free(vals_c3); -	free(errors_c3); +  free(Ls_c); +  free(betas_c); +  free(vals_c1); +  free(errors_c1); +  free(vals_c2); +  free(errors_c2); +  free(vals_c3); +  free(errors_c3); -	return 0; +  return 0;  } - diff --git a/src/corr_test.c b/src/corr_test.c index 01b3e11..cb00212 100644 --- a/src/corr_test.c +++ b/src/corr_test.c @@ -2,26 +2,26 @@  #include "fracture.h"  int main() { -	cholmod_common c; -	CHOL_F(start)(&c); +  cholmod_common c; +  CHOL_F(start)(&c); -	unsigned int width = 64; +  unsigned int width = 64; -	graph_t *network = graph_create(VORONOI_LATTICE, TORUS_BOUND, 128, false); -	net_t *instance = net_create(network, 1e-14, 3, 0, true, &c); -	net_fracture(instance, &c, 1e-10); +  graph_t *network = graph_create(VORONOI_LATTICE, TORUS_BOUND, 128, false); +  net_t *instance = net_create(network, 1e-14, 3, 0, true, &c); +  net_fracture(instance, &c, 1e-10); -	double *corr = get_corr(instance, NULL, &c); +  double *corr = get_corr(instance, NULL, &c); -	for (int i = 0; i < 2 * width; i++) { -		printf("%g ", corr[i]); -	} -	printf("\n"); +  for (int i = 0; i < 2 * width; i++) { +    printf("%g ", corr[i]); +  } +  printf("\n"); -	net_free(instance, &c); -	graph_free(network); +  net_free(instance, &c); +  graph_free(network); -	CHOL_F(finish)(&c); +  CHOL_F(finish)(&c); -	return 0; +  return 0;  } diff --git a/src/fracture.c b/src/fracture.c index 73059ff..ede7a24 100644 --- a/src/fracture.c +++ b/src/fracture.c @@ -3,510 +3,515 @@  #include "fracture.h"  int main(int argc, char *argv[]) { -	int opt; +  int opt; -	// defining variables to be (potentially) set by command line flags -	uint8_t filename_len; -	uint32_t N; -	uint_t L; -	double beta, inf, cutoff, crack_len; -	bool save_data, save_cluster_dist, use_voltage_boundaries, use_dual, save_network, -			 save_crit_stress, save_energy, save_conductivity, -			 save_damage, save_threshold, save_current_load; -	bound_t boundary; -	lattice_t lattice; +  // defining variables to be (potentially) set by command line flags +  uint8_t filename_len; +  uint32_t N; +  uint_t L; +  double beta, inf, cutoff, crack_len; +  bool save_data, save_cluster_dist, use_voltage_boundaries, use_dual, +      save_network, save_crit_stress, save_energy, save_conductivity, +      save_damage, save_threshold, save_current_load; +  bound_t boundary; +  lattice_t lattice; +  // assume filenames will be less than 100 characters -	// assume filenames will be less than 100 characters +  filename_len = 100; -	filename_len = 100; +  // set default values +  N = 100; +  L = 16; +  crack_len = 0.; +  beta = .3; +  inf = 1e10; +  cutoff = 1e-9; +  boundary = FREE_BOUND; +  lattice = VORONOI_LATTICE; +  save_data = false; +  save_cluster_dist = false; +  use_voltage_boundaries = false; +  use_dual = false; +  save_network = false; +  save_crit_stress = false; +  save_damage = false; +  save_conductivity = false; +  save_energy = false; +  save_threshold = false; +  save_current_load = false; -	// set default values +  uint8_t bound_i; +  char boundc2 = 'f'; +  uint8_t lattice_i; +  char lattice_c = 'v'; +  char dual_c = 'o'; -	N = 100; -	L = 16; -	crack_len = 0.; -	beta = .3; -	inf = 1e10; -	cutoff = 1e-9; -	boundary = FREE_BOUND; -	lattice = VORONOI_LATTICE; -	save_data = false; -	save_cluster_dist = false; -	use_voltage_boundaries = false; -	use_dual = false; -	save_network = false; -	save_crit_stress = false; -	save_damage = false; -	save_conductivity = false; -	save_energy = false; -	save_threshold = false; -	save_current_load = false; +  // get commandline options +  while ((opt = getopt(argc, argv, "n:L:b:B:q:dVcoNsCrDl:TE")) != -1) { +    switch (opt) { +    case 'n': +      N = atoi(optarg); +      break; +    case 'L': +      L = atoi(optarg); +      break; +    case 'b': +      beta = atof(optarg); +      break; +    case 'l': +      crack_len = atof(optarg); +      break; +    case 'B': +      bound_i = atoi(optarg); +      switch (bound_i) { +      case 0: +        boundary = FREE_BOUND; +        boundc2 = 'f'; +        break; +      case 1: +        boundary = CYLINDER_BOUND; +        boundc2 = 'c'; +        break; +      case 2: +        boundary = TORUS_BOUND; +        use_voltage_boundaries = true; +        boundc2 = 't'; +        break; +      case 3: +        boundary = EMBEDDED_BOUND; +        boundc2 = 'e'; +        use_dual = true; +        use_voltage_boundaries = true; +        break; +      default: +        printf("boundary specifier must be 0 (FREE_BOUND), 1 (CYLINDER_BOUND), " +               "or 2 (TORUS_BOUND).\n"); +        exit(EXIT_FAILURE); +      } +      break; +    case 'q': +      lattice_i = atoi(optarg); +      switch (lattice_i) { +      case 0: +        lattice = VORONOI_LATTICE; +        lattice_c = 'v'; +        break; +      case 1: +        lattice = DIAGONAL_LATTICE; +        lattice_c = 'd'; +        break; +      case 2: +        lattice = VORONOI_HYPERUNIFORM_LATTICE; +        lattice_c = 'h'; +        break; +      case 3: +        lattice = TRIANGULAR_LATTICE; +        lattice_c = 't'; +        break; +      case 4: +        lattice = SQUARE_LATTICE; +        lattice_c = 's'; +        break; +      default: +        printf("lattice specifier must be 0 (VORONOI_LATTICE), 1 " +               "(DIAGONAL_LATTICE), or 2 (VORONOI_HYPERUNIFORM_LATTICE).\n"); +        exit(EXIT_FAILURE); +      } +      break; +    case 'd': +      save_damage = true; +      break; +    case 'V': +      use_voltage_boundaries = true; +      break; +    case 'D': +      use_dual = true; +      dual_c = 'd'; +      break; +    case 'c': +      save_cluster_dist = true; +      break; +    case 'o': +      save_data = true; +      break; +    case 'N': +      save_network = true; +      break; +    case 's': +      save_crit_stress = true; +      break; +    case 'r': +      save_conductivity = true; +      break; +    case 'E': +      save_energy = true; +      break; +    case 'T': +      save_threshold = true; +      break; +    case 'C': +      save_current_load = true; +      break; +    default: /* '?' */ +      exit(EXIT_FAILURE); +    } +  } -	uint8_t bound_i; -	char boundc2 = 'f'; -	uint8_t lattice_i; -	char lattice_c = 'v'; -	char dual_c = 'o'; +  char boundc; +  if (use_voltage_boundaries) +    boundc = 'v'; +  else +    boundc = 'c'; +  FILE *data_out; +  if (save_data) { +    char *data_filename = (char *)malloc(filename_len * sizeof(char)); +    snprintf(data_filename, filename_len, "data_%c_%c_%c_%c_%u_%g_%g.txt", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +    data_out = fopen(data_filename, "a"); +    free(data_filename); +  } -	// get commandline options +  uint_t max_verts, max_edges; -	while ((opt = getopt(argc, argv, "n:L:b:B:q:dVcoNsCrDl:TE")) != -1) { -		switch (opt) { -			case 'n': -				N = atoi(optarg); -				break; -			case 'L': -				L = atoi(optarg); -				break; -			case 'b': -				beta = atof(optarg); -				break; -			case 'l': -				crack_len = atof(optarg); -				break; -			case 'B': -				bound_i = atoi(optarg); -				switch (bound_i) { -					case 0: -						boundary = FREE_BOUND; -						boundc2 = 'f'; -						break; -					case 1: -						boundary = CYLINDER_BOUND; -						boundc2 = 'c'; -						break; -					case 2: -						boundary = TORUS_BOUND; -						use_voltage_boundaries = true; -						boundc2 = 't'; -						break; -					case 3: -						boundary = EMBEDDED_BOUND; -						boundc2 = 'e'; -						use_dual = true; -						use_voltage_boundaries = true; -						break; -					default: -						printf("boundary specifier must be 0 (FREE_BOUND), 1 (CYLINDER_BOUND), or 2 (TORUS_BOUND).\n"); -						exit(EXIT_FAILURE); -				} -				break; -			case 'q': -				lattice_i = atoi(optarg); -				switch (lattice_i) { -					case 0: -						lattice = VORONOI_LATTICE; -						lattice_c = 'v'; -						break; -					case 1: -						lattice = DIAGONAL_LATTICE; -						lattice_c = 'd'; -						break; -					case 2: -						lattice = VORONOI_HYPERUNIFORM_LATTICE; -						lattice_c = 'h'; -						break; -					case 3: -						lattice = TRIANGULAR_LATTICE; -						lattice_c = 't'; -						break; -          case 4: -            lattice = SQUARE_LATTICE; -            lattice_c = 's'; -            break; -					default: -						printf("lattice specifier must be 0 (VORONOI_LATTICE), 1 (DIAGONAL_LATTICE), or 2 (VORONOI_HYPERUNIFORM_LATTICE).\n"); -						exit(EXIT_FAILURE); -				} -				break; -			case 'd': -				save_damage = true; -				break; -			case 'V': -				use_voltage_boundaries = true; -				break; -			case 'D': -				use_dual = true; -				dual_c = 'd'; -				break; -			case 'c': -				save_cluster_dist = true; -				break; -			case 'o': -				save_data = true; -				break; -			case 'N': -				save_network = true; -				break; -			case 's': -				save_crit_stress = true; -				break; -			case 'r': -				save_conductivity = true; -				break; -			case 'E': -				save_energy = true; -				break; -			case 'T': -				save_threshold = true; -				break; -			case 'C': -				save_current_load = true; -				break; -			default: /* '?' */ -				exit(EXIT_FAILURE); -		} -	} +  // these are very liberal estimates +  max_verts = 4 * pow(L, 2); +  max_edges = 4 * pow(L, 2); +  if (max_verts > CINT_MAX) { +    exit(EXIT_FAILURE); +  } -	char boundc; -	if (use_voltage_boundaries) boundc = 'v'; -	else boundc = 'c'; +  // define arrays for saving cluster and avalanche distributions +  uint32_t *cluster_size_dist; +  uint32_t *avalanche_size_dist; +  char *c_filename; +  char *a_filename; +  if (save_cluster_dist) { +    cluster_size_dist = (uint32_t *)calloc(max_verts, sizeof(uint32_t)); +    avalanche_size_dist = (uint32_t *)calloc(max_edges, sizeof(uint32_t)); -	FILE *data_out; -	if (save_data) { -		char *data_filename = (char *)malloc(filename_len * sizeof(char)); -		snprintf(data_filename, filename_len, "data_%c_%c_%c_%c_%u_%g_%g.txt", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -		data_out = fopen(data_filename, "a"); -		free(data_filename); -	} +    c_filename = (char *)malloc(filename_len * sizeof(char)); +    a_filename = (char *)malloc(filename_len * sizeof(char)); +    snprintf(c_filename, filename_len, "cstr_%c_%c_%c_%c_%d_%g_%g.dat", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +    snprintf(a_filename, filename_len, "avln_%c_%c_%c_%c_%d_%g_%g.dat", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -	uint_t max_verts, max_edges; +    FILE *cluster_out = fopen(c_filename, "rb"); +    FILE *avalanche_out = fopen(a_filename, "rb"); -	// these are very liberal estimates -	max_verts = 4 * pow(L, 2); -	max_edges = 4 * pow(L, 2); +    if (cluster_out != NULL) { +      fread(cluster_size_dist, sizeof(uint32_t), max_verts, cluster_out); +      fclose(cluster_out); +    } +    if (avalanche_out != NULL) { +      fread(avalanche_size_dist, sizeof(uint32_t), max_edges, avalanche_out); +      fclose(avalanche_out); +    } +  } -	if (max_verts > CINT_MAX) { -		exit(EXIT_FAILURE); -	} +  long double *crit_stress; +  if (save_crit_stress) { +    crit_stress = (long double *)malloc(N * sizeof(long double)); +  } -	// define arrays for saving cluster and avalanche distributions -	uint32_t *cluster_size_dist; -	uint32_t *avalanche_size_dist; -	char *c_filename; -	char *a_filename; -	if (save_cluster_dist) { -		cluster_size_dist = -				(uint32_t *)calloc(max_verts, sizeof(uint32_t)); -		avalanche_size_dist = -				(uint32_t *)calloc(max_edges, sizeof(uint32_t)); +  double *conductivity; +  if (save_conductivity) { +    conductivity = (double *)malloc(N * sizeof(double)); +  } -		c_filename = (char *)malloc(filename_len * sizeof(char)); -		a_filename = (char *)malloc(filename_len * sizeof(char)); -		snprintf(c_filename, filename_len, "cstr_%c_%c_%c_%c_%d_%g_%g.dat", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -		snprintf(a_filename, filename_len, "avln_%c_%c_%c_%c_%d_%g_%g.dat", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +  // define arrays for saving damage distributions +  uint32_t *damage; +  char *d_filename; +  if (save_damage) { +    damage = (uint32_t *)calloc(max_edges, sizeof(uint32_t)); -		FILE *cluster_out = fopen(c_filename, "rb"); -		FILE *avalanche_out = fopen(a_filename, "rb"); +    d_filename = (char *)malloc(filename_len * sizeof(char)); +    snprintf(d_filename, filename_len, "damg_%c_%c_%c_%c_%d_%g_%g.dat", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -		if (cluster_out != NULL) { -			fread(cluster_size_dist, sizeof(uint32_t), max_verts, cluster_out); -			fclose(cluster_out); -		} -		if (avalanche_out != NULL) { -			fread(avalanche_size_dist, sizeof(uint32_t), max_edges, avalanche_out); -			fclose(avalanche_out); -		} -	} +    FILE *damage_out = fopen(d_filename, "rb"); -	long double *crit_stress; -	if (save_crit_stress) { -		crit_stress = (long double *)malloc(N * sizeof(long double)); -	} +    if (damage_out != NULL) { +      fread(damage, sizeof(uint32_t), max_edges, damage_out); +      fclose(damage_out); +    } +  } -	double *conductivity; -	if (save_conductivity) { -		conductivity = (double *)malloc(N * sizeof(double)); -	} +  long double *energy; +  if (save_energy) { +    energy = (long double *)malloc(N * sizeof(long double)); +  } -	// define arrays for saving damage distributions -	uint32_t *damage; -	char *d_filename; -	if (save_damage) { -		damage = -				(uint32_t *)calloc(max_edges, sizeof(uint32_t)); +  long double *thresholds; +  if (save_threshold) { +    thresholds = (long double *)malloc(N * sizeof(long double)); +  } -		d_filename = (char *)malloc(filename_len * sizeof(char)); -		snprintf(d_filename, filename_len, "damg_%c_%c_%c_%c_%d_%g_%g.dat", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +  long double *loads; +  if (save_current_load) { +    loads = (long double *)malloc(N * sizeof(long double)); +  } -		FILE *damage_out = fopen(d_filename, "rb"); +  // start cholmod +  cholmod_common c; +  CHOL_F(start)(&c); -		if (damage_out != NULL) { -			fread(damage, sizeof(uint32_t), max_edges, damage_out); -			fclose(damage_out); -		} -	} +  /* if we use voltage boundary conditions, the laplacian matrix is positive +   * definite and we can use a supernodal LL decomposition.  otherwise we need +   * to use the simplicial LDL decomposition +   */ +  if (use_voltage_boundaries) { +    //(&c)->supernodal = CHOLMOD_SUPERNODAL; +    (&c)->supernodal = CHOLMOD_SIMPLICIAL; +  } else { +    (&c)->supernodal = CHOLMOD_SIMPLICIAL; +  } -	long double *energy; -	if (save_energy) { -		energy = (long double *)malloc(N * sizeof(long double)); -	} +  printf("\n"); +  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); -	long double *thresholds; -	if (save_threshold) { -		thresholds = (long double *)malloc(N * sizeof(long double)); -	} +    graph_t *g = graph_create(lattice, boundary, L, use_dual); +    net_t *net = +        net_create(g, inf, beta, crack_len, use_voltage_boundaries, &c); +    net_t *tmp_net = net_copy(net, &c); +    data_t *data = net_fracture(tmp_net, &c, cutoff); +    net_free(tmp_net, &c); -	long double *loads; -	if (save_current_load) { -		loads = (long double *)malloc(N * sizeof(long double)); -	} +    uint_t max_pos = 0; +    long double max_val = 0; +    double cond0; +    { +      double *tmp_voltages = net_voltages(net, &c); +      cond0 = net_conductivity(net, tmp_voltages); +      free(tmp_voltages); +    } -	// start cholmod -	cholmod_common c; -	CHOL_F(start)(&c); +    for (uint_t j = 0; j < data->num_broken; j++) { +      long double val = data->extern_field[j]; -	/* if we use voltage boundary conditions, the laplacian matrix is positive -	 * definite and we can use a supernodal LL decomposition.  otherwise we need -	 * to use the simplicial LDL decomposition -	 */ -	if (use_voltage_boundaries) { -		//(&c)->supernodal = CHOLMOD_SUPERNODAL; -		(&c)->supernodal = CHOLMOD_SIMPLICIAL; -	} else { -		(&c)->supernodal = CHOLMOD_SIMPLICIAL; -	} +      if (val > max_val) { +        max_pos = j; +        max_val = val; +      } +    } +    uint_t av_size = 0; +    long double cur_val = 0; -	printf("\n"); -	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); +    for (uint_t j = 0; j < max_pos; j++) { +      uint_t next_broken = data->break_list[j]; -		graph_t *g = graph_create(lattice, boundary, L, use_dual); -		net_t *net = net_create(g, inf, beta, crack_len, use_voltage_boundaries, &c); -		net_t *tmp_net = net_copy(net, &c); -		data_t *data = net_fracture(tmp_net, &c, cutoff); -		net_free(tmp_net, &c); +      break_edge(net, next_broken, &c); -		uint_t max_pos = 0; -		long double max_val = 0; +      long double val = data->extern_field[j]; +      if (save_cluster_dist) { +        if (val < cur_val) { +          av_size++; +        } -		double cond0; -		{ -			double *tmp_voltages = net_voltages(net, &c); -			cond0 = net_conductivity(net, tmp_voltages); -			free(tmp_voltages); -		} +        if (val > cur_val) { +          avalanche_size_dist[av_size]++; +          av_size = 0; +          cur_val = val; +        } +      } +    } -		for (uint_t j = 0; j < data->num_broken; j++) { -			long double val = data->extern_field[j]; +    if (save_crit_stress) +      crit_stress[i] = data->extern_field[max_pos]; -			if (val > max_val) { -				max_pos = j; -				max_val = val; -			} -		} +    if (save_conductivity) { +      if (max_pos > 0) { +        conductivity[i] = data->conductivity[max_pos - 1]; +      } else { +        conductivity[i] = cond0; +      } +    } -		uint_t av_size = 0; -		long double cur_val = 0; +    if (save_damage) { +      uint_t would_break = 0; +      double *tmp_voltage = net_voltages(net, &c); +      double *tmp_current = net_currents(net, tmp_voltage, &c); +      free(tmp_voltage); +      for (uint_t j = 0; j < g->ne; j++) { +        bool broken = net->fuses[j]; +        bool under_thres = +            net->thres[j] < net->thres[data->break_list[max_pos]]; +        bool zero_field = fabs(tmp_current[j]) < cutoff; +        if (!broken && under_thres && zero_field) { +          break_edge(net, j, &c); +        } +      } +      damage[net->num_broken]++; +      free(tmp_current); +    } -		for (uint_t j = 0; j < max_pos; j++) { -			uint_t next_broken = data->break_list[j]; +    if (save_energy) { +      long double tmp_energy = 0; +      if (max_pos > 0) { +        long double sigma1 = data->extern_field[0]; +        double cond1 = cond0; +        for (uint_t j = 0; j < max_pos - 1; j++) { +          long double sigma2 = data->extern_field[j + 1]; +          double cond2 = data->conductivity[j]; +          if (sigma2 > sigma1) { +            tmp_energy += 0.5 * gsl_pow_2(sigma1) * (1 - cond2 / cond1) / cond1; +            sigma1 = sigma2; +            cond1 = cond2; +          } +        } +      } +      energy[i] = tmp_energy; +    } -			break_edge(net, next_broken, &c); +    if (save_threshold) { +      thresholds[i] = net->thres[data->break_list[max_pos]]; +    } -			long double val = data->extern_field[j]; -			if (save_cluster_dist) { -				if (val < cur_val) { -					av_size++; -				} +    if (save_current_load) { +      loads[i] = +          data->extern_field[max_pos] / net->thres[data->break_list[max_pos]]; +    } -				if (val > cur_val) { -					avalanche_size_dist[av_size]++; -					av_size = 0; -					cur_val = val; -				} -			} -		} +    if (save_data) { +      for (uint_t j = 0; j < data->num_broken; j++) { +        fprintf(data_out, "%u %Lg %g ", data->break_list[j], +                data->extern_field[j], data->conductivity[j]); +      } +      fprintf(data_out, "\n"); +    } -		if (save_crit_stress) crit_stress[i] = data->extern_field[max_pos]; +    data_free(data); +    if (save_network) { +      FILE *net_out = fopen("network.txt", "w"); +      for (uint_t j = 0; j < g->nv; j++) { +        fprintf(net_out, "%f %f ", g->vx[2 * j], g->vx[2 * j + 1]); +      } +      fprintf(net_out, "\n"); +      for (uint_t j = 0; j < g->ne; j++) { +        fprintf(net_out, "%u %u ", g->ev[2 * j], g->ev[2 * j + 1]); +      } +      fprintf(net_out, "\n"); +      for (uint_t j = 0; j < g->dnv; j++) { +        fprintf(net_out, "%f %f ", g->dvx[2 * j], g->dvx[2 * j + 1]); +      } +      fprintf(net_out, "\n"); +      for (uint_t j = 0; j < g->ne; j++) { +        fprintf(net_out, "%u %u ", g->dev[2 * j], g->dev[2 * j + 1]); +      } +      fprintf(net_out, "\n"); +      for (uint_t j = 0; j < g->ne; j++) { +        fprintf(net_out, "%d ", net->fuses[j]); +      } +      fclose(net_out); +    } -		if (save_conductivity) { -			if (max_pos > 0) { -				conductivity[i] = data->conductivity[max_pos - 1]; -			} else { -				conductivity[i] = cond0; -			} -		} +    if (save_cluster_dist) { +      uint_t *tmp_cluster_dist = get_cluster_dist(net); +      for (uint_t j = 0; j < g->dnv; j++) { +        cluster_size_dist[j] += tmp_cluster_dist[j]; +      } +      free(tmp_cluster_dist); +    } -		if (save_damage) { -			uint_t would_break = 0; -			double *tmp_voltage = net_voltages(net, &c); -			double *tmp_current = net_currents(net, tmp_voltage, &c); -			free(tmp_voltage); -			for (uint_t j = 0; j < g->ne; j++) { -				bool broken = net->fuses[j]; -				bool under_thres = net->thres[j] < net->thres[data->break_list[max_pos]]; -				bool zero_field = fabs(tmp_current[j]) < cutoff; -				if (!broken && under_thres && zero_field) { -					break_edge(net, j, &c); -				} -			} -			damage[net->num_broken]++; -			free(tmp_current); -		} +    net_free(net, &c); +    graph_free(g); +  } -		if (save_energy) { -			long double tmp_energy = 0; -			if (max_pos > 0) { -				long double sigma1 = data->extern_field[0]; -				double cond1 = cond0; -				for (uint_t j = 0; j < max_pos - 1; j++) { -					long double sigma2 = data->extern_field[j+1]; -					double cond2 = data->conductivity[j]; -					if (sigma2 > sigma1) { -						tmp_energy += 0.5 * gsl_pow_2(sigma1) * (1 - cond2 / cond1) / cond1; -						sigma1 = sigma2; -						cond1 = cond2; -					} -				} -			} -			energy[i] = tmp_energy; -		} +  printf("\033[F\033[JFRACTURE: COMPLETE\n"); -		if (save_threshold) { -			thresholds[i] = net->thres[data->break_list[max_pos]]; -		} +  if (save_cluster_dist) { +    FILE *cluster_out = fopen(c_filename, "wb"); +    FILE *avalanche_out = fopen(a_filename, "wb"); -		if (save_current_load) { -			loads[i] = data->extern_field[max_pos] / net->thres[data->break_list[max_pos]]; -		} +    fwrite(cluster_size_dist, sizeof(uint32_t), max_verts, cluster_out); +    fwrite(avalanche_size_dist, sizeof(uint32_t), max_edges, avalanche_out); -		if (save_data) { -			for (uint_t j = 0; j < data->num_broken; j++) { -				fprintf(data_out, "%u %Lg %g ", data->break_list[j], -								data->extern_field[j], data->conductivity[j]); -			} -			fprintf(data_out, "\n"); -		} +    fclose(cluster_out); +    fclose(avalanche_out); -		data_free(data); -		if (save_network) { -			FILE *net_out = fopen("network.txt", "w"); -			for (uint_t j = 0; j < g->nv; j++) { -				fprintf(net_out, "%f %f ", g->vx[2 * j], -						g->vx[2 * j + 1]); -			} -			fprintf(net_out, "\n"); -			for (uint_t j = 0; j < g->ne; j++) { -				fprintf(net_out, "%u %u ", g->ev[2 * j], g->ev[2 * j + 1]); -			} -			fprintf(net_out, "\n"); -			for (uint_t j = 0; j < g->dnv; j++) { -				fprintf(net_out, "%f %f ", g->dvx[2 * j], -						g->dvx[2 * j + 1]); -			} -			fprintf(net_out, "\n"); -			for (uint_t j = 0; j < g->ne; j++) { -				fprintf(net_out, "%u %u ", g->dev[2 * j], g->dev[2 * j + 1]); -			} -			fprintf(net_out, "\n"); -			for (uint_t j = 0; j < g->ne; j++) { -				fprintf(net_out, "%d ", net->fuses[j]); -			} -			fclose(net_out); -		} +    free(c_filename); +    free(a_filename); +    free(cluster_size_dist); +    free(avalanche_size_dist); +  } -		if (save_cluster_dist) { -			uint_t *tmp_cluster_dist = get_cluster_dist(net); -			for (uint_t j = 0; j < g->dnv; j++) { -				cluster_size_dist[j] += tmp_cluster_dist[j]; -			} -			free(tmp_cluster_dist); -		} +  if (save_conductivity) { +    char *cond_filename = (char *)malloc(filename_len * sizeof(char)); +    snprintf(cond_filename, filename_len, "cond_%c_%c_%c_%c_%d_%g_%g.dat", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +    FILE *cond_file = fopen(cond_filename, "ab"); +    fwrite(conductivity, sizeof(double), N, cond_file); +    fclose(cond_file); +    free(cond_filename); +    free(conductivity); +  } +  if (save_energy) { +    char *tough_filename = (char *)malloc(filename_len * sizeof(char)); +    snprintf(tough_filename, filename_len, "enrg_%c_%c_%c_%c_%d_%g_%g.dat", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +    FILE *tough_file = fopen(tough_filename, "ab"); +    fwrite(energy, sizeof(long double), N, tough_file); +    fclose(tough_file); +    free(tough_filename); +    free(energy); +  } -		net_free(net, &c); -		graph_free(g); -	} +  if (save_threshold) { +    char *thres_filename = (char *)malloc(filename_len * sizeof(char)); +    snprintf(thres_filename, filename_len, "thrs_%c_%c_%c_%c_%d_%g_%g.dat", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +    FILE *thres_file = fopen(thres_filename, "ab"); +    fwrite(thresholds, sizeof(long double), N, thres_file); +    fclose(thres_file); +    free(thres_filename); +    free(thresholds); +  } -	printf("\033[F\033[JFRACTURE: COMPLETE\n"); +  if (save_current_load) { +    char *load_filename = (char *)malloc(filename_len * sizeof(char)); +    snprintf(load_filename, filename_len, "load_%c_%c_%c_%c_%d_%g_%g.dat", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +    FILE *load_file = fopen(load_filename, "ab"); +    fwrite(loads, sizeof(long double), N, load_file); +    fclose(load_file); +    free(load_filename); +    free(loads); +  } -	if (save_cluster_dist) { -		FILE *cluster_out = fopen(c_filename, "wb"); -		FILE *avalanche_out = fopen(a_filename, "wb"); +  if (save_damage) { +    FILE *hdam_file = fopen(d_filename, "wb"); +    fwrite(damage, sizeof(uint32_t), max_edges, hdam_file); +    fclose(hdam_file); +    free(d_filename); +    free(damage); +  } -		fwrite(cluster_size_dist, sizeof(uint32_t), max_verts, cluster_out); -		fwrite(avalanche_size_dist, sizeof(uint32_t), max_edges, avalanche_out); +  if (save_data) { +    fclose(data_out); +  } -		fclose(cluster_out); -		fclose(avalanche_out); +  if (save_crit_stress) { +    char *str_filename = (char *)malloc(filename_len * sizeof(char)); +    snprintf(str_filename, filename_len, "strs_%c_%c_%c_%c_%d_%g_%g.dat", +             lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); +    FILE *str_file = fopen(str_filename, "ab"); +    fwrite(crit_stress, sizeof(long double), N, str_file); +    fclose(str_file); +    free(str_filename); +    free(crit_stress); +  } -		free(c_filename); -		free(a_filename); -		free(cluster_size_dist); -		free(avalanche_size_dist); -	} +  CHOL_F(finish)(&c); -	if (save_conductivity) { -		char *cond_filename = (char *)malloc(filename_len * sizeof(char)); -		snprintf(cond_filename, filename_len, "cond_%c_%c_%c_%c_%d_%g_%g.dat", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -		FILE *cond_file = fopen(cond_filename, "ab"); -		fwrite(conductivity, sizeof(double), N, cond_file); -		fclose(cond_file); -		free(cond_filename); -		free(conductivity); -	} - -	if (save_energy) { -		char *tough_filename = (char *)malloc(filename_len * sizeof(char)); -		snprintf(tough_filename, filename_len, "enrg_%c_%c_%c_%c_%d_%g_%g.dat", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -		FILE *tough_file = fopen(tough_filename, "ab"); -		fwrite(energy, sizeof(long double), N, tough_file); -		fclose(tough_file); -		free(tough_filename); -		free(energy); -	} - -	if (save_threshold) { -		char *thres_filename = (char *)malloc(filename_len * sizeof(char)); -		snprintf(thres_filename, filename_len, "thrs_%c_%c_%c_%c_%d_%g_%g.dat", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -		FILE *thres_file = fopen(thres_filename, "ab"); -		fwrite(thresholds, sizeof(long double), N, thres_file); -		fclose(thres_file); -		free(thres_filename); -		free(thresholds); -	} - -	if (save_current_load) { -		char *load_filename = (char *)malloc(filename_len * sizeof(char)); -		snprintf(load_filename, filename_len, "load_%c_%c_%c_%c_%d_%g_%g.dat", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -		FILE *load_file = fopen(load_filename, "ab"); -		fwrite(loads, sizeof(long double), N, load_file); -		fclose(load_file); -		free(load_filename); -		free(loads); -	} - -	if (save_damage) { -		FILE *hdam_file = fopen(d_filename, "wb"); -		fwrite(damage, sizeof(uint32_t), max_edges, hdam_file); -		fclose(hdam_file); -		free(d_filename); -		free(damage); -	} - -	if (save_data) { -		fclose(data_out); -	} - -	if (save_crit_stress) { -		char *str_filename = (char *)malloc(filename_len * sizeof(char)); -		snprintf(str_filename, filename_len, "strs_%c_%c_%c_%c_%d_%g_%g.dat", lattice_c, dual_c, boundc, boundc2, L, beta, crack_len); -		FILE *str_file = fopen(str_filename, "ab"); -		fwrite(crit_stress, sizeof(long double), N, str_file); -		fclose(str_file); -		free(str_filename); -		free(crit_stress); -	} - -	CHOL_F(finish)(&c); - -	return 0; +  return 0;  } diff --git a/src/fracture.h b/src/fracture.h index b1114fb..7b4f146 100644 --- a/src/fracture.h +++ b/src/fracture.h @@ -17,7 +17,6 @@  #include <stdbool.h>  #include <stdint.h>  #include <stdio.h> -#include <stdio.h>  #include <stdlib.h>  #include <string.h>  #include <sys/types.h> @@ -26,7 +25,6 @@  #include <jst/graph.h>  #include <jst/rand.h> -  // these defs allow me to switch to long int cholmod in a sitch  #define int_t int  #define uint_t unsigned int @@ -36,68 +34,69 @@  #define GSL_RAND_GEN gsl_rng_mt19937  typedef struct { -	const graph_t *graph; -	bool *fuses; -	long double *thres; -	double inf; -	cholmod_dense *boundary_cond; -	cholmod_factor *factor; -	bool voltage_bound; -	uint_t num_broken; -	uint_t dim; -	uint_t nep; -	uint_t *evp; -	cholmod_sparse *voltcurmat; +  const graph_t *graph; +  bool *fuses; +  long double *thres; +  double inf; +  cholmod_dense *boundary_cond; +  cholmod_factor *factor; +  bool voltage_bound; +  uint_t num_broken; +  uint_t dim; +  uint_t nep; +  uint_t *evp; +  cholmod_sparse *voltcurmat;  } net_t;  typedef struct { -	uint_t num_broken; -	uint_t *break_list; -	double *conductivity; -	long double *extern_field; +  uint_t num_broken; +  uint_t *break_list; +  double *conductivity; +  long double *extern_field;  } data_t; -intptr_t *run_voronoi(uint_t num_coords, double *coords, bool periodic, double xmin, double xmax, double ymin, double ymax); +intptr_t *run_voronoi(uint_t num_coords, double *coords, bool periodic, +                      double xmin, double xmax, double ymin, double ymax); -cholmod_sparse *gen_adjacency(const net_t *net, bool dual, bool use_gp, bool symmetric, cholmod_common *c); +cholmod_sparse *gen_adjacency(const net_t *net, bool dual, bool use_gp, +                              bool symmetric, cholmod_common *c);  cholmod_sparse *gen_laplacian(const net_t *net, cholmod_common *c); -int edge_to_verts(uint_t width, bool periodic, uint_t edge, -									bool index); +int edge_to_verts(uint_t width, bool periodic, uint_t edge, bool index); -int dual_edge_to_verts(uint_t width, bool periodic, uint_t edge, -											 bool index); +int dual_edge_to_verts(uint_t width, bool periodic, uint_t edge, bool index); -double dual_vert_to_coord(uint_t width, bool periodic, uint_t vert, -													bool index); +double dual_vert_to_coord(uint_t width, bool periodic, uint_t vert, bool index); -void factor_update(cholmod_factor *factor, uint_t v1, uint_t v2, cholmod_common *c); +void factor_update(cholmod_factor *factor, uint_t v1, uint_t v2, +                   cholmod_common *c);  void factor_update2(cholmod_factor *factor, uint_t v1, cholmod_common *c);  void net_notch(net_t *net, double notch_len, cholmod_common *c);  data_t *net_fracture(net_t *net, cholmod_common *c, double cutoff);  double *net_voltages(const net_t *net, cholmod_common *c); -double *net_currents(const net_t *net, const double *voltages, cholmod_common *c); +double *net_currents(const net_t *net, const double *voltages, +                     cholmod_common *c);  double net_conductivity(const net_t *net, const double *voltages);  void update_boundary(net_t *instance, const double *avg_field); -FILE *get_file(const char *prefix, uint_t width, uint_t crack, -							 double beta, uint_t iter, uint_t num_iter, -							 uint_t num, bool read); +FILE *get_file(const char *prefix, uint_t width, uint_t crack, double beta, +               uint_t iter, uint_t num_iter, uint_t num, bool read);  double update_beta(double beta, uint_t width, const double *stress, -									 const double *damage, double bound_total); +                   const double *damage, double bound_total);  cholmod_sparse *gen_voltcurmat(uint_t num_edges, uint_t num_verts, -															 uint_t *edges_to_verts, cholmod_common *c); +                               uint_t *edges_to_verts, cholmod_common *c);  net_t *net_copy(const net_t *net, cholmod_common *c);  void net_free(net_t *instance, cholmod_common *c); -net_t *net_create(const graph_t *g, double inf, double beta, double notch_len, bool vb, cholmod_common *c); +net_t *net_create(const graph_t *g, double inf, double beta, double notch_len, +                  bool vb, cholmod_common *c);  bool break_edge(net_t *instance, uint_t edge, cholmod_common *c); @@ -112,11 +111,13 @@ double *get_corr(net_t *instance, uint_t **dists, cholmod_common *c);  double *bin_values(graph_t *network, uint_t width, double *values); -cholmod_dense *bound_set(const graph_t *g, bool vb, double notch_len, cholmod_common *c); +cholmod_dense *bound_set(const graph_t *g, bool vb, double notch_len, +                         cholmod_common *c);  data_t *data_create(uint_t num_edges);  void data_free(data_t *data); -void data_update(data_t *data, uint_t last_broke, long double strength, double conductivity); +void data_update(data_t *data, uint_t last_broke, long double strength, +                 double conductivity);  long double rand_dist_pow(const gsl_rng *r, double beta); diff --git a/src/long_anal_process.c b/src/long_anal_process.c index ba29152..d4ec4e0 100644 --- a/src/long_anal_process.c +++ b/src/long_anal_process.c @@ -1,156 +1,158 @@  #include "fracture.h" +#include <gsl/gsl_blas.h> +#include <gsl/gsl_matrix.h>  #include <gsl/gsl_sf_erf.h>  #include <gsl/gsl_sf_laguerre.h> -#include <gsl/gsl_matrix.h>  #include <gsl/gsl_vector.h> -#include <gsl/gsl_blas.h>  #include <sys/stat.h>  void get_mean(uint_t len, long double *data, long double result[2]) { -	long double total = 0; +  long double total = 0; -	for (uint_t i = 0; i < len; i++) { -		total += data[i]; -	} +  for (uint_t i = 0; i < len; i++) { +    total += data[i]; +  } -	long double mean = total / len; -	long double total_sq_diff = 0; +  long double mean = total / len; +  long double total_sq_diff = 0; -	for (uint_t i = 0; i < len; i++) { -		total_sq_diff += pow(mean - data[i], 2); -	} +  for (uint_t i = 0; i < len; i++) { +    total_sq_diff += pow(mean - data[i], 2); +  } -	long double mean_err = sqrt(total_sq_diff) / len; +  long double mean_err = sqrt(total_sq_diff) / len; -	result[0] = mean; -	result[1] = mean_err; +  result[0] = mean; +  result[1] = mean_err;  } -void get_mom(uint_t len, uint_t n, long double *data, long double mean[2], long double result[2]) { -	long double total_n_diff = 0; -	long double total_n2_diff = 0; +void get_mom(uint_t len, uint_t n, long double *data, long double mean[2], +             long double result[2]) { +  long double total_n_diff = 0; +  long double total_n2_diff = 0; -	for (uint_t i = 0; i < len; i++) { -		total_n_diff += pow(fabsl(mean[0] - data[i]), n); -		total_n2_diff += pow(fabsl(mean[0] - data[i]), 2 * n); -	} +  for (uint_t i = 0; i < len; i++) { +    total_n_diff += pow(fabsl(mean[0] - data[i]), n); +    total_n2_diff += pow(fabsl(mean[0] - data[i]), 2 * n); +  } -	long double mom = total_n_diff / len; -	long double mom_err = sqrt(total_n2_diff) / len; +  long double mom = total_n_diff / len; +  long double mom_err = sqrt(total_n2_diff) / len; -	result[0] = mom; -	result[1] = mom_err; +  result[0] = mom; +  result[1] = mom_err;  }  int main(int argc, char *argv[]) { -	uint_t nc = argc - 1; -	uint_t *Ls_c = (uint_t *)malloc(nc * sizeof(uint_t)); -	double *betas_c = (double *)malloc(nc * sizeof(double)); -	long double *vals_c1 = (long double *)malloc(nc * sizeof(long double)); -	long double *errors_c1 = (long double *)malloc(nc * sizeof(long double)); -	long double *vals_c2 = (long double *)malloc(nc * sizeof(long double)); -	long double *errors_c2 = (long double *)malloc(nc * sizeof(long double)); -	long double *vals_c3 = (long double *)malloc(nc * sizeof(long double)); -	long double *errors_c3 = (long double *)malloc(nc * sizeof(long double)); -	long double *vals_c4 = (long double *)malloc(nc * sizeof(long double)); -	long double *errors_c4 = (long double *)malloc(nc * sizeof(long double)); - -	char *out_filename = (char *)malloc(100 * sizeof(char)); -	uint_t out_filename_len = 0; +  uint_t nc = argc - 1; +  uint_t *Ls_c = (uint_t *)malloc(nc * sizeof(uint_t)); +  double *betas_c = (double *)malloc(nc * sizeof(double)); +  long double *vals_c1 = (long double *)malloc(nc * sizeof(long double)); +  long double *errors_c1 = (long double *)malloc(nc * sizeof(long double)); +  long double *vals_c2 = (long double *)malloc(nc * sizeof(long double)); +  long double *errors_c2 = (long double *)malloc(nc * sizeof(long double)); +  long double *vals_c3 = (long double *)malloc(nc * sizeof(long double)); +  long double *errors_c3 = (long double *)malloc(nc * sizeof(long double)); +  long double *vals_c4 = (long double *)malloc(nc * sizeof(long double)); +  long double *errors_c4 = (long double *)malloc(nc * sizeof(long double)); -	for (uint_t i = 0; i < nc; i++) { -		char *fn = argv[1 + i]; -		char *buff = (char *)malloc(20 * sizeof(char)); -		uint_t pos = 0; uint_t c = 0; -		while (c < 5) { -			if (fn[pos] == '_') { -				c++; -			} -			if (i == 0) { -				out_filename[pos] = fn[pos]; -				out_filename_len++; -			} -			pos++; -		} -		c = 0; -		while (fn[pos] != '_') { -			buff[c] = fn[pos]; -			pos++; -			c++; -		} -		buff[c] = '\0'; -		Ls_c[i] = atoi(buff); -		c = 0; -		pos++; -		while (fn[pos] != '_') { -			buff[c] = fn[pos]; -			pos++; -			c++; -		} -		buff[c] = '\0'; -		betas_c[i] = atof(buff); -		free(buff); +  char *out_filename = (char *)malloc(100 * sizeof(char)); +  uint_t out_filename_len = 0; -		struct stat info; -		stat(fn, &info); -		uint_t num_bytes = info.st_size; -		uint_t dist_len = (num_bytes * sizeof(char)) / sizeof(long double); +  for (uint_t i = 0; i < nc; i++) { +    char *fn = argv[1 + i]; +    char *buff = (char *)malloc(20 * sizeof(char)); +    uint_t pos = 0; +    uint_t c = 0; +    while (c < 5) { +      if (fn[pos] == '_') { +        c++; +      } +      if (i == 0) { +        out_filename[pos] = fn[pos]; +        out_filename_len++; +      } +      pos++; +    } +    c = 0; +    while (fn[pos] != '_') { +      buff[c] = fn[pos]; +      pos++; +      c++; +    } +    buff[c] = '\0'; +    Ls_c[i] = atoi(buff); +    c = 0; +    pos++; +    while (fn[pos] != '_') { +      buff[c] = fn[pos]; +      pos++; +      c++; +    } +    buff[c] = '\0'; +    betas_c[i] = atof(buff); +    free(buff); -		long double *dist = malloc(dist_len * sizeof(long double)); -		FILE *dist_file = fopen(fn, "rb"); -		fread(dist, sizeof(long double), dist_len, dist_file); -		fclose(dist_file); -		{ -			long double mom1[2]; -			get_mean(dist_len, dist, mom1); -			vals_c1[i] = mom1[0]; -			errors_c1[i] = mom1[1]; +    struct stat info; +    stat(fn, &info); +    uint_t num_bytes = info.st_size; +    uint_t dist_len = (num_bytes * sizeof(char)) / sizeof(long double); -			long double mom2[2]; -			get_mom(dist_len, 2, dist, mom1, mom2); -			vals_c2[i] = mom2[0]; -			errors_c2[i] = mom2[1]; +    long double *dist = malloc(dist_len * sizeof(long double)); +    FILE *dist_file = fopen(fn, "rb"); +    fread(dist, sizeof(long double), dist_len, dist_file); +    fclose(dist_file); +    { +      long double mom1[2]; +      get_mean(dist_len, dist, mom1); +      vals_c1[i] = mom1[0]; +      errors_c1[i] = mom1[1]; -			long double mom3[2]; -			get_mom(dist_len, 3, dist, mom1, mom3); -			vals_c3[i] = mom3[0]; -			errors_c3[i] = mom3[1]; +      long double mom2[2]; +      get_mom(dist_len, 2, dist, mom1, mom2); +      vals_c2[i] = mom2[0]; +      errors_c2[i] = mom2[1]; -			long double mom4[2]; -			get_mom(dist_len, 4, dist, mom1, mom4); -			vals_c4[i] = mom4[0]; -			errors_c4[i] = mom4[1]; -		} -		free(dist); -	} +      long double mom3[2]; +      get_mom(dist_len, 3, dist, mom1, mom3); +      vals_c3[i] = mom3[0]; +      errors_c3[i] = mom3[1]; -	out_filename[out_filename_len-1] = '.'; -	out_filename[out_filename_len] = 't'; -	out_filename[out_filename_len+1] = 'x'; -	out_filename[out_filename_len+2] = 't'; -	out_filename[out_filename_len+3] = '\0'; +      long double mom4[2]; +      get_mom(dist_len, 4, dist, mom1, mom4); +      vals_c4[i] = mom4[0]; +      errors_c4[i] = mom4[1]; +    } +    free(dist); +  } -	FILE *out_file = fopen(out_filename, "w"); +  out_filename[out_filename_len - 1] = '.'; +  out_filename[out_filename_len] = 't'; +  out_filename[out_filename_len + 1] = 'x'; +  out_filename[out_filename_len + 2] = 't'; +  out_filename[out_filename_len + 3] = '\0'; -	for (uint_t i = 0; i < nc; i++) { -		fprintf(out_file, "%u %g %Lg %Lg %Lg %Lg %Lg %Lg %Lg %Lg\n", Ls_c[i], betas_c[i], vals_c1[i], errors_c1[i], vals_c2[i], errors_c2[i], vals_c3[i], errors_c3[i], vals_c4[i], errors_c4[i]); -	} +  FILE *out_file = fopen(out_filename, "w"); -	fclose(out_file); +  for (uint_t i = 0; i < nc; i++) { +    fprintf(out_file, "%u %g %Lg %Lg %Lg %Lg %Lg %Lg %Lg %Lg\n", Ls_c[i], +            betas_c[i], vals_c1[i], errors_c1[i], vals_c2[i], errors_c2[i], +            vals_c3[i], errors_c3[i], vals_c4[i], errors_c4[i]); +  } +  fclose(out_file); -	free(Ls_c); -	free(betas_c); -	free(vals_c1); -	free(errors_c1); -	free(vals_c2); -	free(errors_c2); -	free(vals_c3); -	free(errors_c3); -	free(vals_c4); -	free(errors_c4); +  free(Ls_c); +  free(betas_c); +  free(vals_c1); +  free(errors_c1); +  free(vals_c2); +  free(errors_c2); +  free(vals_c3); +  free(errors_c3); +  free(vals_c4); +  free(errors_c4); -	return 0; +  return 0;  } -  | 
