parallelized lensing kernel
@@ -274,9 +274,6 @@ void ccl_get_lensing_mag_kernel(ccl_cosmology *cosmo,
int *status) {
ccl_f1d_t *nz_f = NULL;
ccl_f1d_t *sz_f = NULL;
- gsl_integration_workspace *w;
- integ_lensing_pars *ipar;
- int local_status;
// Prepare N(z) spline
nz_f = ccl_f1d_t_new(nz, z_arr, nz_arr, 0, 0,
@@ -313,14 +310,19 @@ void ccl_get_lensing_mag_kernel(ccl_cosmology *cosmo,
}
}
- double lens_prefac = get_lensing_prefactor(cosmo, status);
- double chi, a, z, mgfac;
- int ichi;
+ if(*status==0) {
+ #pragma omp parallel default(none) \
+ shared(cosmo, z_max, i_nz_norm, sz_f, nz_f, \
+ nchi, chi_arr, wL_arr, status)
+ {
+ double chi, a, z, mgfac, lens_prefac;
+ int ichi, local_status;
+ integ_lensing_pars *ipar = NULL;
+ gsl_integration_workspace *w = NULL;
- w = NULL;
- ipar = NULL;
local_status = *status;
+ lens_prefac = get_lensing_prefactor(cosmo, &local_status);
if (local_status == 0) {
ipar = malloc(sizeof(integ_lensing_pars));
w = gsl_integration_workspace_alloc(cosmo->gsl_params.N_ITERATION);
@@ -340,6 +342,7 @@ void ccl_get_lensing_mag_kernel(ccl_cosmology *cosmo,
}
//Populate arrays
+ #pragma omp for
for (ichi=0; ichi < nchi; ichi++) {
if (local_status == 0) {
chi = chi_arr[ichi];
@@ -356,14 +359,17 @@ void ccl_get_lensing_mag_kernel(ccl_cosmology *cosmo,
} else {
wL_arr[ichi] = NAN;
}
- }
+ } //end omp for
gsl_integration_workspace_free(w);
free(ipar);
- if (local_status != 0) {
+ if (local_status) {
+ #pragma omp atomic write
*status = CCL_ERROR_INTEG;
}
+ } //end omp parallel
+ }
ccl_f1d_t_free(nz_f);
ccl_f1d_t_free(sz_f);
