| #include "XSbench_header.h" |
| |
| #ifdef MPI |
| #include<mpi.h> |
| #endif |
| |
| // Generates randomized energy grid for each nuclide |
| // Note that this is done as part of initialization (serial), so |
| // rand() is used. |
| void generate_grids( NuclideGridPoint ** nuclide_grids, |
| long n_isotopes, long n_gridpoints ) { |
| for( long i = 0; i < n_isotopes; i++ ) |
| for( long j = 0; j < n_gridpoints; j++ ) |
| { |
| nuclide_grids[i][j].energy =((double)rand()/(double)RAND_MAX); |
| nuclide_grids[i][j].total_xs =((double)rand()/(double)RAND_MAX); |
| nuclide_grids[i][j].elastic_xs =((double)rand()/(double)RAND_MAX); |
| nuclide_grids[i][j].absorbtion_xs=((double)rand()/(double)RAND_MAX); |
| nuclide_grids[i][j].fission_xs =((double)rand()/(double)RAND_MAX); |
| nuclide_grids[i][j].nu_fission_xs=((double)rand()/(double)RAND_MAX); |
| } |
| } |
| |
| // Verification version of this function (tighter control over RNG) |
| void generate_grids_v( NuclideGridPoint ** nuclide_grids, |
| long n_isotopes, long n_gridpoints ) { |
| for( long i = 0; i < n_isotopes; i++ ) |
| for( long j = 0; j < n_gridpoints; j++ ) |
| { |
| nuclide_grids[i][j].energy = rn_v(); |
| nuclide_grids[i][j].total_xs = rn_v(); |
| nuclide_grids[i][j].elastic_xs = rn_v(); |
| nuclide_grids[i][j].absorbtion_xs= rn_v(); |
| nuclide_grids[i][j].fission_xs = rn_v(); |
| nuclide_grids[i][j].nu_fission_xs= rn_v(); |
| } |
| } |
| |
| // Sorts the nuclide grids by energy (lowest -> highest) |
| void sort_nuclide_grids( NuclideGridPoint ** nuclide_grids, long n_isotopes, |
| long n_gridpoints ) |
| { |
| int (*cmp) (const void *, const void *); |
| cmp = NGP_compare; |
| |
| for( long i = 0; i < n_isotopes; i++ ) |
| qsort( nuclide_grids[i], n_gridpoints, sizeof(NuclideGridPoint), |
| cmp ); |
| |
| // error debug check |
| /* |
| for( int i = 0; i < n_isotopes; i++ ) |
| { |
| printf("NUCLIDE %d ==============================\n", i); |
| for( int j = 0; j < n_gridpoints; j++ ) |
| printf("E%d = %lf\n", j, nuclide_grids[i][j].energy); |
| } |
| */ |
| } |
| |
| // Allocates unionized energy grid, and assigns union of energy levels |
| // from nuclide grids to it. |
| GridPoint * generate_energy_grid( long n_isotopes, long n_gridpoints, |
| NuclideGridPoint ** nuclide_grids) { |
| int mype = 0; |
| |
| #ifdef MPI |
| MPI_Comm_rank(MPI_COMM_WORLD, &mype); |
| #endif |
| |
| if( mype == 0 ) printf("Generating Unionized Energy Grid...\n"); |
| |
| long n_unionized_grid_points = n_isotopes*n_gridpoints; |
| int (*cmp) (const void *, const void *); |
| cmp = NGP_compare; |
| |
| GridPoint * energy_grid = (GridPoint *)malloc( n_unionized_grid_points |
| * sizeof( GridPoint ) ); |
| if( mype == 0 ) printf("Copying and Sorting all nuclide grids...\n"); |
| |
| NuclideGridPoint ** n_grid_sorted = gpmatrix( n_isotopes, n_gridpoints ); |
| |
| |
| memcpy( n_grid_sorted[0], nuclide_grids[0], n_isotopes*n_gridpoints* |
| sizeof( NuclideGridPoint ) ); |
| |
| qsort( &n_grid_sorted[0][0], n_unionized_grid_points, |
| sizeof(NuclideGridPoint), cmp); |
| |
| if( mype == 0 ) printf("Assigning energies to unionized grid...\n"); |
| |
| for( long i = 0; i < n_unionized_grid_points; i++ ) |
| energy_grid[i].energy = n_grid_sorted[0][i].energy; |
| |
| |
| gpmatrix_free(n_grid_sorted); |
| |
| int * full = (int *) malloc( n_isotopes * n_unionized_grid_points |
| * sizeof(int) ); |
| if( full == NULL ) |
| { |
| fprintf(stderr,"ERROR - Out Of Memory!\n"); |
| exit(1); |
| } |
| |
| for( long i = 0; i < n_unionized_grid_points; i++ ) |
| energy_grid[i].xs_ptrs = &full[n_isotopes * i]; |
| |
| // debug error checking |
| /* |
| for( int i = 0; i < n_unionized_grid_points; i++ ) |
| printf("E%d = %lf\n", i, energy_grid[i].energy); |
| */ |
| |
| return energy_grid; |
| } |
| |
| // Searches each nuclide grid for the closest energy level and assigns |
| // pointer from unionized grid to the correct spot in the nuclide grid. |
| // This process is time consuming, as the number of binary searches |
| // required is: binary searches = n_gridpoints * n_isotopes^2 |
| void set_grid_ptrs( GridPoint * energy_grid, NuclideGridPoint ** nuclide_grids, |
| long n_isotopes, long n_gridpoints ) |
| { |
| int mype = 0; |
| |
| #ifdef MPI |
| MPI_Comm_rank(MPI_COMM_WORLD, &mype); |
| #endif |
| |
| if( mype == 0 ) printf("Assigning pointers to Unionized Energy Grid...\n"); |
| #ifdef OPENMP |
| #pragma omp parallel for default(none) \ |
| shared( energy_grid, nuclide_grids, n_isotopes, n_gridpoints, mype ) |
| #endif |
| for( long i = 0; i < n_isotopes * n_gridpoints ; i++ ) |
| { |
| int nthreads = 1, tid = 0; |
| double quarry = energy_grid[i].energy; |
| |
| #ifdef OPENMP |
| nthreads = omp_get_num_threads(); |
| tid = omp_get_thread_num(); |
| #endif |
| |
| if( INFO && mype == 0 && tid == 0 && i % 200 == 0 ) |
| printf("\rAligning Unionized Grid...(%.0lf%% complete)", |
| 100.0 * (double) i / (n_isotopes*n_gridpoints / |
| nthreads) ); |
| for( long j = 0; j < n_isotopes; j++ ) |
| { |
| // j is the nuclide i.d. |
| // log n binary search |
| energy_grid[i].xs_ptrs[j] = |
| binary_search( nuclide_grids[j], quarry, n_gridpoints); |
| } |
| } |
| if( mype == 0 ) printf("\n"); |
| |
| //test |
| /* |
| for( int i=0; i < n_isotopes * n_gridpoints; i++ ) |
| for( int j = 0; j < n_isotopes; j++ ) |
| printf("E = %.4lf\tNuclide %d->%p->%.4lf\n", |
| energy_grid[i].energy, |
| j, |
| energy_grid[i].xs_ptrs[j], |
| (energy_grid[i].xs_ptrs[j])->energy |
| ); |
| */ |
| } |