48 #ifdef MEASURE_TIME_FFTW 60 static void flags_cp(NFFT(plan) *dst, NFFT(plan) *src)
63 dst->f_hat = src->f_hat;
67 dst->my_fftw_plan1 = src->my_fftw_plan1;
68 dst->my_fftw_plan2 = src->my_fftw_plan2;
71 static void time_accuracy(
int d,
int N,
int M,
int n,
int m,
unsigned test_ndft,
72 unsigned test_pre_full_psi)
80 NFFT(plan) p_pre_phi_hut;
82 NFFT(plan) p_pre_lin_psi;
83 NFFT(plan) p_pre_fg_psi;
85 NFFT(plan) p_pre_full_psi;
87 printf("%d\t%d\t", d, N);
89 for (r = 0; r < d; r++)
97 swapndft = (C*) NFFT(malloc)((size_t)(M) *
sizeof(C));
99 NFFT(init_guru)(&p, d, NN, M, nn, m,
100 MALLOC_X | MALLOC_F_HAT | MALLOC_F |
101 FFTW_INIT | FFT_OUT_OF_PLACE,
102 FFTW_MEASURE | FFTW_DESTROY_INPUT);
105 NFFT(vrand_shifted_unit_double)(p.x, p.d * p.M_total);
107 NFFT(init_guru)(&p_pre_phi_hut, d, NN, M, nn, m, PRE_PHI_HUT, 0);
108 flags_cp(&p_pre_phi_hut, &p);
109 NFFT(precompute_one_psi)(&p_pre_phi_hut);
113 NFFT(init_guru)(&p_fg_psi, d, NN, M, nn, m, FG_PSI, 0);
114 flags_cp(&p_fg_psi, &p);
115 NFFT(precompute_one_psi)(&p_fg_psi);
118 NFFT(init_guru)(&p_pre_lin_psi, d, NN, M, nn, m, PRE_LIN_PSI, 0);
119 flags_cp(&p_pre_lin_psi, &p);
120 NFFT(precompute_one_psi)(&p_pre_lin_psi);
124 NFFT(init_guru)(&p_pre_fg_psi, d, NN, M, nn, m, PRE_FG_PSI, 0);
125 flags_cp(&p_pre_fg_psi, &p);
126 NFFT(precompute_one_psi)(&p_pre_fg_psi);
129 NFFT(init_guru)(&p_pre_psi, d, NN, M, nn, m, PRE_PSI, 0);
130 flags_cp(&p_pre_psi, &p);
131 NFFT(precompute_one_psi)(&p_pre_psi);
133 if (test_pre_full_psi)
135 NFFT(init_guru)(&p_pre_full_psi, d, NN, M, nn, m, PRE_FULL_PSI, 0);
136 flags_cp(&p_pre_full_psi, &p);
137 NFFT(precompute_one_psi)(&p_pre_full_psi);
141 NFFT(vrand_unit_complex)(p.f_hat, p.N_total);
146 CSWAP(p.f, swapndft);
150 while (t_ndft < K(0.01))
154 NFFT(trafo_direct)(&p);
156 t = NFFT(elapsed_seconds)(t1, t0);
161 CSWAP(p.f, swapndft);
168 NFFT(trafo)(&p_pre_phi_hut);
170 NFFT(trafo)(&p_fg_psi);
172 p_fg_psi.MEASURE_TIME_t[2] = MKNAN(
"");
173 NFFT(trafo)(&p_pre_lin_psi);
175 NFFT(trafo)(&p_pre_fg_psi);
177 p_pre_fg_psi.MEASURE_TIME_t[2] = MKNAN(
"");
178 NFFT(trafo)(&p_pre_psi);
179 if (test_pre_full_psi)
180 NFFT(trafo)(&p_pre_full_psi);
182 p_pre_full_psi.MEASURE_TIME_t[2] = MKNAN(
"");
185 p.MEASURE_TIME_t[1] = MKNAN(
"");
188 e = NFFT(error_l_2_complex)(swapndft, p.f, p.M_total);
193 "%.2" __FES__
"\t%d\t%.2" __FES__
"\t%.2" __FES__
"\t%.2" __FES__
"\t%.2" __FES__
"\t%.2" __FES__
"\t%.2" __FES__
"\t%.2" __FES__
"\t%.2" __FES__
"\t%.2" __FES__
"\t%.2" __FES__
"\n",
194 t_ndft, m, e, p.MEASURE_TIME_t[0], p_pre_phi_hut.MEASURE_TIME_t[0],
195 p.MEASURE_TIME_t[1], p.MEASURE_TIME_t[2], p_fg_psi.MEASURE_TIME_t[2],
196 p_pre_lin_psi.MEASURE_TIME_t[2], p_pre_fg_psi.MEASURE_TIME_t[2],
197 p_pre_psi.MEASURE_TIME_t[2], p_pre_full_psi.MEASURE_TIME_t[2]);
202 if (test_pre_full_psi)
203 NFFT(finalize)(&p_pre_full_psi);
204 NFFT(finalize)(&p_pre_psi);
206 NFFT(finalize)(&p_pre_fg_psi);
207 NFFT(finalize)(&p_pre_lin_psi);
209 NFFT(finalize)(&p_fg_psi);
210 NFFT(finalize)(&p_pre_phi_hut);
214 NFFT(free)(swapndft);
217 static void accuracy_pre_lin_psi(
int d,
int N,
int M,
int n,
int m,
int K)
225 for (r = 0; r < d; r++)
232 swapndft = (C*) NFFT(malloc)((size_t)(M) *
sizeof(C));
234 NFFT(init_guru)(&p, d, NN, M, nn, m,
235 MALLOC_X | MALLOC_F_HAT | MALLOC_F |
236 PRE_PHI_HUT | PRE_LIN_PSI |
237 FFTW_INIT | FFT_OUT_OF_PLACE,
238 FFTW_MEASURE | FFTW_DESTROY_INPUT);
243 p.psi = (R*) NFFT(malloc)((size_t)((p.K + 1) * p.d) *
sizeof(R));
246 NFFT(precompute_one_psi)(&p);
249 NFFT(vrand_shifted_unit_double)(p.x, p.d * p.M_total);
252 NFFT(vrand_unit_complex)(p.f_hat, p.N_total);
255 CSWAP(p.f, swapndft);
256 NFFT(trafo_direct)(&p);
257 CSWAP(p.f, swapndft);
261 e = NFFT(error_l_2_complex)(swapndft, p.f, p.M_total);
264 printf(
"$%.1" __FES__
"$&\t", e);
270 NFFT(free)(swapndft);
273 int main(
int argc,
char **argv)
279 fprintf(stderr,
"flags type first last trials d m\n");
283 if ((test == 0) && (atoi(argv[1]) < 2))
285 fprintf(stderr,
"MEASURE_TIME in infft.h not set\n");
289 fprintf(stderr,
"Testing different precomputation schemes for the nfft.\n");
290 fprintf(stderr,
"Columns: d, N=M, t_ndft, e_nfft, t_D, t_pre_phi_hut, ");
291 fprintf(stderr,
"t_fftw, t_B, t_fg_psi, t_pre_lin_psi, t_pre_fg_psi, ");
292 fprintf(stderr,
"t_pre_psi, t_pre_full_psi\n\n");
294 int arg2 = atoi(argv[2]);
295 int arg3 = atoi(argv[3]);
296 int arg4 = atoi(argv[4]);
299 if (atoi(argv[1]) == 0)
301 int d = atoi(argv[5]);
302 int m = atoi(argv[6]);
304 for (l = arg2; l <= arg3; l++)
306 int N = (int)(1U << l);
307 int M = (int)(1U << (d * l));
308 for (trial = 0; trial < arg4; trial++)
310 time_accuracy(d, N, M, 2 * N, m, 0, 0);
314 else if (atoi(argv[1]) == 1)
316 int d = atoi(argv[5]);
317 int N = atoi(argv[6]);
320 for (m = arg2; m <= arg3; m++)
322 for (trial = 0; trial < arg4; trial++)
324 time_accuracy(d, N, (
int)(LRINT(POW((R)(N), (R)(d)))), 2 * N, m, 1, 1);
328 else if (atoi(argv[1]) == 2)
330 int d = atoi(argv[5]);
331 int N = atoi(argv[6]);
332 int m = atoi(argv[7]);
334 printf(
"$\\log_2(K/(m+1))$&\t");
336 for (l = arg2; l < arg3; l++)
337 printf(
"$%d$&\t", l);
339 printf(
"$%d$\\\\\n", arg3);
341 printf(
"$\\tilde E_2$&\t");
342 for (l = arg2; l <= arg3; l++)
344 int x = (m + 1) * (
int)(1U << l);
345 accuracy_pre_lin_psi(d, N, (
int)(LRINT(POW((R)(N), (R)(d)))), 2 * N, m, x);
#define CSWAP(x, y)
Swap two vectors.