diff --git a/phigrape.cpp b/phigrape.cpp index a82b749..1b89885 100644 --- a/phigrape.cpp +++ b/phigrape.cpp @@ -314,46 +314,6 @@ double L[3]; // needed in pn_bh_spin.c // double var_sort[N_MAX]; -double3 xcm, vcm, mom, - xdc, vdc, - a2, a3, a2dot1; - -char processor_name[MPI_MAX_PROCESSOR_NAME], - inp_fname[30], - out_fname[30]; - -/* global variables */ - -int N, N_star, N_bh, - ind[N_MAX], name[N_MAX]; -double m[N_MAX], pot[N_MAX], t[N_MAX], dt[N_MAX]; -double3 x[N_MAX], v[N_MAX], a[N_MAX], adot[N_MAX]; - -/* local variables */ - -int n_loc, ind_loc[N_MAX_loc]; -double m_loc[N_MAX_loc], pot_loc[N_MAX_loc], t_loc[N_MAX_loc], dt_loc[N_MAX_loc]; -double3 x_loc[N_MAX_loc], v_loc[N_MAX_loc], - a_loc[N_MAX_loc], adot_loc[N_MAX_loc]; - -/* data for active particles */ - -int n_act, ind_act[N_MAX]; -double m_act[N_MAX], - pot_act[N_MAX], t_act[N_MAX], dt_act[N_MAX], - pot_act_new[N_MAX], - pot_act_tmp[N_MAX], - pot_act_tmp_loc[N_MAX]; - -double3 x_act[N_MAX], v_act[N_MAX], - a_act[N_MAX], adot_act[N_MAX], - x_act_new[N_MAX], v_act_new[N_MAX], - a_act_new[N_MAX], adot_act_new[N_MAX], - a_act_tmp[N_MAX], adot_act_tmp[N_MAX], - a_act_tmp_loc[N_MAX], adot_act_tmp_loc[N_MAX]; - -FILE *inp, *out; - double CPU_time_real0, CPU_time_user0, CPU_time_syst0; double CPU_time_real, CPU_time_user, CPU_time_syst; @@ -584,12 +544,11 @@ gettimeofday(&tv, NULL); } -void read_data(int *diskstep, int *N, double *time_cur, int ind[], double m[], double3 x[], double3 v[]) +void read_data(char inp_fname[], int *diskstep, int *N, double *time_cur, int ind[], double m[], double3 x[], double3 v[]) { auto inp = fopen(inp_fname, "r"); fscanf(inp, "%d \n", diskstep); fscanf(inp, "%d \n", N); - printf("zzzzzz reading *%s* %d\n", inp_fname, *N); fscanf(inp, "%lE \n", time_cur); for (int i=0; i<*N; i++) fscanf(inp,"%d %lE %lE %lE %lE %lE %lE %lE \n", &ind[i], &m[i], &x[i][0], &x[i][1], &x[i][2], &v[i][0], &v[i][1], &v[i][2]); fclose(inp); @@ -597,7 +556,7 @@ void read_data(int *diskstep, int *N, double *time_cur, int ind[], double m[], d void write_snap_data(char out_fname[], int diskstep, int N, double time_cur, int ind[], double m[], double3 x[], double3 v[]) { - out = fopen(out_fname, "w"); + auto out = fopen(out_fname, "w"); fprintf(out,"%06d \n", diskstep); fprintf(out,"%07d \n", N); fprintf(out,"%.10E \n", time_cur); @@ -611,28 +570,7 @@ void write_snap_data(char out_fname[], int diskstep, int N, double time_cur, int fclose(out); } -// // // // // // // // // // // void write_cont_data() // virtually identical to write_snap_data() -// // // // // // // // // // // { -// // // // // // // // // // // out = fopen("data.con","w"); -// // // // // // // // // // // fprintf(out,"%06d \n", diskstep); -// // // // // // // // // // // -// // // // // // // // // // // fprintf(out,"%07d \n", N); -// // // // // // // // // // // -// // // // // // // // // // // fprintf(out,"%.16E \n", time_cur); -// // // // // // // // // // // -// // // // // // // // // // // for (i=0; ilive_smbh_count == 2) { @@ -735,7 +673,7 @@ void write_bh_data(double time_cur, double m[], double3 x[], double3 v[]) fprintf(out,"\n"); fclose(out); } else if (config->live_smbh_count == 1) { - out = fopen("bh.dat", "a"); + auto out = fopen("bh.dat", "a"); double tmp_r = sqrt( SQR(x[0][0]) + SQR(x[0][1]) + SQR(x[0][2]) ); double tmp_v = sqrt( SQR(v[0][0]) + SQR(v[0][1]) + SQR(v[0][2]) ); double tmp_a = sqrt( SQR(a[0][0]) + SQR(a[0][1]) + SQR(a[0][2]) ); @@ -753,13 +691,13 @@ void write_bh_data(double time_cur, double m[], double3 x[], double3 v[]) } } -void write_bh_nb_data(double time_cur, double m[], double3 x[], double3 v[]) +void write_bh_nb_data(double time_cur, int N, double m[], double3 x[], double3 v[]) { int nb = config->live_smbh_neighbor_number; int ind_sort[N_MAX]; double var_sort[N_MAX]; - out = fopen("bh_neighbors.dat", "a"); + auto out = fopen("bh_neighbors.dat", "a"); /* 1st BH */ @@ -1039,7 +977,13 @@ for (i=0; ilive_smbh_output) write_bh_data(time_cur, m, x, v); + if (config->live_smbh_output) write_bh_data(time_cur, m, x, v, pot, a, adot, dt); /* Write BH NB data... */ - if (config->live_smbh_neighbor_output) write_bh_nb_data(time_cur, m, x, v); + if (config->live_smbh_neighbor_output) write_bh_nb_data(time_cur, N, m, x, v); } /* if (myRank == rootRank) */ @@ -2261,7 +2253,13 @@ double dt_disk, dt_contr, t_disk=0.0, t_contr=0.0, DT_ACT_PRED += (CPU_tmp_user - CPU_tmp_user0); #endif - calc_self_grav(min_t, eps2, g6_calls); + calc_self_grav(min_t, eps2, g6_calls, n_loc, + n_act, ind_act, + x_act_new, v_act_new, + pot_act_tmp, + a_act_tmp, + adot_act_tmp, + h2_i); /* Reduce the "global" vectors from "local" on all the nodes */ @@ -2580,10 +2578,10 @@ double dt_disk, dt_contr, t_disk=0.0, t_contr=0.0, if (time_cur >= t_bh) { if (myRank == rootRank) { /* Write BH data... */ - if (config->live_smbh_output) write_bh_data(time_cur, m, x, v); + if (config->live_smbh_output) write_bh_data(time_cur, m, x, v, pot, a, adot, dt); /* Write BH NB data... */ - if (config->live_smbh_neighbor_output) write_bh_nb_data(time_cur, m, x, v); + if (config->live_smbh_neighbor_output) write_bh_nb_data(time_cur, N, m, x, v); } /* if (myRank == rootRank) */ @@ -2593,7 +2591,7 @@ double dt_disk, dt_contr, t_disk=0.0, t_contr=0.0, if (time_cur >= t_contr) { if (myRank == rootRank) { - energy_contr(time_cur, timesteps, n_act_sum, g6_calls, rcm_sum, vcm_sum, E_tot_0, E_tot_corr_0, E_tot_corr_sd_0, skip_con); + energy_contr(time_cur, timesteps, n_act_sum, g6_calls, rcm_sum, vcm_sum, E_tot_0, E_tot_corr_0, E_tot_corr_sd_0, skip_con, N, m, x, v, pot); /* write cont data */