33 void bench_openmp_readfile(FILE *infile,
int *trafo_adjoint,
int *N,
int *M,
double **x, C **f_hat, C **f)
39 fscanf(infile,
"%d %d %d", trafo_adjoint, N, M);
41 *f_hat = (C*)
nfft_malloc((2*(*N)+2) * (2*(*N)+2) *
sizeof(C));
44 memset(*f_hat,0U,(2*(*N)+2) * (2*(*N)+2) *
sizeof(C));
45 memset(*f,0U,(*M)*
sizeof(C));
48 fftw_import_wisdom_from_filename(
"nfsft_benchomp_detail_threads.plan");
50 fftw_import_wisdom_from_filename(
"nfsft_benchomp_detail_single.plan");
53 nfsft_init_guru(&plan, *N, *M, NFSFT_MALLOC_X | NFSFT_MALLOC_F |
54 NFSFT_MALLOC_F_HAT | NFSFT_NORMALIZED | NFSFT_PRESERVE_F_HAT,
55 PRE_PHI_HUT | FFTW_INIT | FFT_OUT_OF_PLACE, 6);
58 fftw_export_wisdom_to_filename(
"nfsft_benchomp_detail_threads.plan");
60 fftw_export_wisdom_to_filename(
"nfsft_benchomp_detail_single.plan");
63 for (j=0; j < *M; j++)
65 fscanf(infile,
"%lg", (*x)+2*j+t);
69 for (k = 0; k <= *N; k++)
70 for (n = -k; n <= k; n++)
72 fscanf(infile,
"%lg %lg", &re, &im);
73 (*f_hat)[NFSFT_INDEX(k,n,&plan)] = re + _Complex_I * im;
78 for (j=0; j < *M; j++)
80 fscanf(infile,
"%lg %lg", &re, &im);
81 (*f)[j] = re + _Complex_I * im;
85 nfsft_finalize(&plan);
88 void bench_openmp(
int trafo_adjoint,
int N,
int M,
double *x, C *f_hat, C *f,
int m,
int nfsft_flags,
int psi_flags)
95 double tt_total, tt_pre;
114 nfsft_init_guru(&plan, N, M, nfsft_flags | NFSFT_MALLOC_X | NFSFT_MALLOC_F |
115 NFSFT_MALLOC_F_HAT | NFSFT_NORMALIZED | NFSFT_PRESERVE_F_HAT,
116 PRE_PHI_HUT | psi_flags | FFTW_INIT | FFT_OUT_OF_PLACE, m);
124 for (j=0; j < plan.
M_total; j++)
126 for (t=0; t < 2; t++)
128 plan.
x[2*j+t] = x[2*j+t];
131 if (trafo_adjoint==0)
133 memset(plan.
f_hat,0U,plan.
N_total*
sizeof(
double _Complex));
134 for (k = 0; k <= plan.
N; k++)
135 for (n = -k; n <= k; n++)
139 plan.
f_hat[NFSFT_INDEX(k,n,&plan)] = f_hat[NFSFT_INDEX(k,n,&plan)];
144 for (j=0; j < plan.
M_total; j++)
151 memset(plan.
f_hat,0U,plan.
N_total*
sizeof(
double _Complex));
156 nfsft_precompute_x(&plan);
160 if (trafo_adjoint==0)
163 nfsft_adjoint(&plan);
172 #ifndef MEASURE_TIME_FFTW 179 nfsft_finalize(&plan);
182 int main(
int argc,
char **argv)
184 int m, nfsft_flags, psi_flags;
186 int trafo_adjoint, N, M, r;
195 nthreads = atoi(argv[5]);
197 omp_set_num_threads(nthreads);
204 nfsft_flags = atoi(argv[2]);
205 psi_flags = atoi(argv[3]);
206 nrepeat = atoi(argv[4]);
208 bench_openmp_readfile(stdin, &trafo_adjoint, &N, &M, &x, &f_hat, &f);
211 nfsft_precompute(N,1000.0,0U,0U);
213 for (r = 0; r < nrepeat; r++)
214 bench_openmp(trafo_adjoint, N, M, x, f_hat, f, m, nfsft_flags, psi_flags);
double * x
the nodes for ,
fftw_complex * f_hat
Fourier coefficients.
R nfft_elapsed_seconds(ticks t1, ticks t0)
Return number of elapsed seconds between two time points.
data structure for an NFSFT (nonequispaced fast spherical Fourier transform) plan with double precisi...
NFFT_INT M_total
Total number of samples.
double MEASURE_TIME_t[3]
Measured time for each step if MEASURE_TIME is set.
void * nfft_malloc(size_t n)
NFFT_INT N_total
Total number of Fourier coefficients.