00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #include <stdio.h>
00022 #include <math.h>
00023 #include <string.h>
00024 #include <stdlib.h>
00025 #include <complex.h>
00026
00027 #include "nfft3util.h"
00028 #include "nfft3.h"
00029
00030 int global_n;
00031 int global_d;
00032
00033 static int comp1(const void *x,const void *y)
00034 {
00035 return ((* (double*) x)<(* (double*) y)?-1:1);
00036 }
00037
00038 static int comp2(const void *x,const void *y)
00039 {
00040 int nx0,nx1,ny0,ny1;
00041 nx0=global_n*(* ((double*)x+0));
00042 nx1=global_n*(* ((double*)x+1));
00043 ny0=global_n*(* ((double*)y+0));
00044 ny1=global_n*(* ((double*)y+1));
00045
00046 if(nx0<ny0)
00047 return -1;
00048 else
00049 if(nx0==ny0)
00050 if(nx1<ny1)
00051 return -1;
00052 else
00053 return 1;
00054 else
00055 return 1;
00056 }
00057
00058 static int comp3(const void *x,const void *y)
00059 {
00060 int nx0,nx1,nx2,ny0,ny1,ny2;
00061 nx0=global_n*(* ((double*)x+0));
00062 nx1=global_n*(* ((double*)x+1));
00063 nx2=global_n*(* ((double*)x+2));
00064 ny0=global_n*(* ((double*)y+0));
00065 ny1=global_n*(* ((double*)y+1));
00066 ny2=global_n*(* ((double*)y+2));
00067
00068 if(nx0<ny0)
00069 return -1;
00070 else
00071 if(nx0==ny0)
00072 if(nx1<ny1)
00073 return -1;
00074 else
00075 if(nx1==ny1)
00076 if(nx2<ny2)
00077 return -1;
00078 else
00079 return 1;
00080 else
00081 return 1;
00082 else
00083 return 1;
00084 }
00085
00086 void measure_time_nfft(int d, int N, unsigned test_ndft)
00087 {
00088 int r, M, NN[d], nn[d],j;
00089 double t, t_fft, t_ndft, t_nfft;
00090
00091 nfft_plan p;
00092 fftw_plan p_fft;
00093
00094 double auxC=pow(2,29);
00095
00096 printf("\\verb+%d+&\t",(int)(log(N)/log(2)*d+0.5));
00097
00098 for(r=0,M=1;r<d;r++)
00099 {
00100 M=N*M;
00101 NN[r]=N;
00102 nn[r]=2*N;
00103 }
00104
00105 nfft_init_guru(&p, d, NN, M, nn, 2,
00106 PRE_PHI_HUT|
00107 PRE_FULL_PSI| MALLOC_F_HAT| MALLOC_X| MALLOC_F|
00108 FFTW_INIT| FFT_OUT_OF_PLACE,
00109 FFTW_MEASURE| FFTW_DESTROY_INPUT);
00110
00111 p_fft=fftw_plan_dft(d, NN, p.f_hat, p.f, FFTW_FORWARD, FFTW_MEASURE);
00112
00114 nfft_vrand_shifted_unit_double(p.x, p.d*p.M_total);
00115
00116 global_n=nn[0];
00117 global_d=d;
00118 switch(d)
00119 {
00120 case 1: { qsort(p.x,p.M_total,d*sizeof(double),comp1); break; }
00121 case 2: { qsort(p.x,p.M_total,d*sizeof(double),comp2); break; }
00122 case 3: { qsort(p.x,p.M_total,d*sizeof(double),comp3); break; }
00123 }
00124
00125 nfft_precompute_one_psi(&p);
00126
00128 nfft_vrand_unit_complex(p.f_hat, p.N_total);
00129
00131 t_fft=0;
00132 r=0;
00133 while(t_fft<1)
00134 {
00135 r++;
00136 t=nfft_second();
00137 fftw_execute(p_fft);
00138 t=nfft_second()-t;
00139 t_fft+=t;
00140 }
00141 t_fft/=r;
00142
00143
00144 printf("\\verb+%.1e+ &\t",t_fft);
00145
00147 if(test_ndft)
00148 {
00149 t_ndft=0;
00150 r=0;
00151 while(t_ndft<1)
00152 {
00153 r++;
00154 t=nfft_second();
00155 ndft_trafo(&p);
00156 t=nfft_second()-t;
00157 t_ndft+=t;
00158 }
00159 t_ndft/=r;
00160
00161 printf("\\verb+%.1e+ &\t",t_ndft);
00162 }
00163 else
00164
00165 printf("\\verb+*+\t&\t");
00166
00168 t_nfft=0;
00169 r=0;
00170 while(t_nfft<1)
00171 {
00172 r++;
00173 t=nfft_second();
00174 switch(d)
00175 {
00176 case 1: { nfft_trafo_1d(&p); break; }
00177 case 2: { nfft_trafo_2d(&p); break; }
00178 case 3: { nfft_trafo_3d(&p); break; }
00179 default: nfft_trafo(&p);
00180 }
00181 t=nfft_second()-t;
00182 t_nfft+=t;
00183 }
00184 t_nfft/=r;
00185
00186
00187 printf("\\verb+%.1e+ & \\verb+(%3.1f)+\\\\\n",t_nfft,t_nfft/t_fft);
00188
00189 fftw_destroy_plan(p_fft);
00190 nfft_finalize(&p);
00191 }
00192
00193 void measure_time_nfft_XXX2(int d, int N, unsigned test_ndft)
00194 {
00195 int r, M, NN[d], nn[d];
00196 double t, t_fft, t_ndft, t_nfft;
00197
00198 nfft_plan p;
00199 fftw_plan p_fft;
00200
00201 printf("%d\t",(int)(log(N)/log(2)*d+0.5)); fflush(stdout);
00202
00203 for(r=0,M=1;r<d;r++)
00204 {
00205 M=N*M;
00206 NN[r]=N;
00207 nn[r]=2*N;
00208 }
00209
00210 nfft_init_guru(&p, d, NN, M, nn, 2,
00211 PRE_PHI_HUT|
00212 PRE_PSI|
00213 MALLOC_F_HAT| MALLOC_X| MALLOC_F|
00214 FFTW_INIT| FFT_OUT_OF_PLACE,
00215 FFTW_MEASURE| FFTW_DESTROY_INPUT);
00216
00217 p_fft=fftw_plan_dft(d, NN, p.f_hat, p.f, FFTW_FORWARD, FFTW_MEASURE);
00218
00219 double _Complex *swapndft=(double _Complex*)nfft_malloc(p.M_total*sizeof(double _Complex));
00220
00222 nfft_vrand_shifted_unit_double(p.x, p.d*p.M_total);
00223
00224 qsort(p.x,p.M_total,d*sizeof(double),comp1);
00225
00226
00227 nfft_precompute_one_psi(&p);
00228
00230 nfft_vrand_unit_complex(p.f_hat, p.N_total);
00231
00233 t_fft=0;
00234 r=0;
00235 while(t_fft<0.1)
00236 {
00237 r++;
00238 t=nfft_second();
00239 fftw_execute(p_fft);
00240 t=nfft_second()-t;
00241 t_fft+=t;
00242 }
00243 t_fft/=r;
00244
00245 printf("%.1e\t",t_fft);
00246
00248 if(test_ndft)
00249 {
00250 NFFT_SWAP_complex(p.f,swapndft);
00251 t_ndft=0;
00252 r=0;
00253 while(t_ndft<0.1)
00254 {
00255 r++;
00256 t=nfft_second();
00257 ndft_trafo(&p);
00258 t=nfft_second()-t;
00259 t_ndft+=t;
00260 }
00261 t_ndft/=r;
00262 printf("%.1e\t",t_ndft);
00263 NFFT_SWAP_complex(p.f,swapndft);
00264 }
00265 else
00266 printf("\t");
00267
00269 t_nfft=0;
00270 r=0;
00271 while(t_nfft<0.1)
00272 {
00273 r++;
00274 t=nfft_second();
00275 nfft_trafo(&p);
00276 t=nfft_second()-t;
00277 t_nfft+=t;
00278 }
00279 t_nfft/=r;
00280 printf("%.1e\t",t_nfft);
00281 if(test_ndft)
00282 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f, p.M_total));
00283
00285 t_nfft=0;
00286 r=0;
00287 while(t_nfft<0.1)
00288 {
00289 r++;
00290 t=nfft_second();
00291 nfft_trafo_1d(&p);
00292 t=nfft_second()-t;
00293 t_nfft+=t;
00294 }
00295 t_nfft/=r;
00296 printf("%.1e\t",t_nfft);
00297 if(test_ndft)
00298 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f, p.M_total));
00299
00300 printf("\n");
00301
00302 nfft_free(swapndft);
00303 fftw_destroy_plan(p_fft);
00304 nfft_finalize(&p);
00305 }
00306
00307 void measure_time_nfft_XXX3(int d, int N, unsigned test_ndft)
00308 {
00309 int r, M, NN[d], nn[d];
00310 double t, t_fft, t_ndft, t_nfft;
00311
00312 nfft_plan p;
00313 fftw_plan p_fft;
00314
00315 printf("%d\t",(int)(log(N)/log(2)*d+0.5)); fflush(stdout);
00316
00317 for(r=0,M=1;r<d;r++)
00318 {
00319 M=N*M;
00320 NN[r]=N;
00321 nn[r]=2*N;
00322 }
00323
00324 nfft_init_guru(&p, d, NN, M, nn, 2,
00325 PRE_PHI_HUT|
00326 PRE_PSI|
00327 MALLOC_F_HAT| MALLOC_X| MALLOC_F|
00328 FFTW_INIT| FFT_OUT_OF_PLACE,
00329 FFTW_MEASURE| FFTW_DESTROY_INPUT);
00330
00331 p_fft=fftw_plan_dft(d, NN, p.f, p.f_hat, FFTW_BACKWARD, FFTW_MEASURE);
00332
00333 double _Complex *swapndft=(double _Complex*)nfft_malloc(p.N_total*sizeof(double _Complex));
00334
00336 nfft_vrand_shifted_unit_double(p.x, p.d*p.M_total);
00337
00338 qsort(p.x,p.M_total,d*sizeof(double),comp1);
00339
00340
00341 nfft_precompute_one_psi(&p);
00342
00344 nfft_vrand_unit_complex(p.f, p.N_total);
00345
00347 t_fft=0;
00348 r=0;
00349 while(t_fft<0.1)
00350 {
00351 r++;
00352 t=nfft_second();
00353 fftw_execute(p_fft);
00354 t=nfft_second()-t;
00355 t_fft+=t;
00356 }
00357 t_fft/=r;
00358
00359 printf("%.1e\t",t_fft);
00360
00362 if(test_ndft)
00363 {
00364 NFFT_SWAP_complex(p.f_hat,swapndft);
00365 t_ndft=0;
00366 r=0;
00367 while(t_ndft<0.1)
00368 {
00369 r++;
00370 t=nfft_second();
00371 ndft_adjoint(&p);
00372 t=nfft_second()-t;
00373 t_ndft+=t;
00374 }
00375 t_ndft/=r;
00376 printf("%.1e\t",t_ndft);
00377 NFFT_SWAP_complex(p.f_hat,swapndft);
00378 }
00379 else
00380 printf("\t");
00381
00383 t_nfft=0;
00384 r=0;
00385 while(t_nfft<0.1)
00386 {
00387 r++;
00388 t=nfft_second();
00389 nfft_adjoint(&p);
00390 t=nfft_second()-t;
00391 t_nfft+=t;
00392 }
00393 t_nfft/=r;
00394 printf("%.1e\t",t_nfft);
00395 if(test_ndft)
00396 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f_hat, p.N_total));
00397
00399 t_nfft=0;
00400 r=0;
00401 while(t_nfft<0.1)
00402 {
00403 r++;
00404 t=nfft_second();
00405 nfft_adjoint_1d(&p);
00406 t=nfft_second()-t;
00407 t_nfft+=t;
00408 }
00409 t_nfft/=r;
00410 printf("%.1e\t",t_nfft);
00411 if(test_ndft)
00412 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f_hat, p.N_total));
00413
00414 printf("\n");
00415
00416 nfft_free(swapndft);
00417 fftw_destroy_plan(p_fft);
00418 nfft_finalize(&p);
00419 }
00420
00421
00422
00423
00424 void measure_time_nfft_XXX4(int d, int N, unsigned test_ndft)
00425 {
00426 int j,r, M, NN[d], nn[d];
00427 double t, t_fft, t_ndft, t_nfft;
00428
00429 nfft_plan p;
00430 fftw_plan p_fft;
00431
00432 printf("%d\t",(int)(log(N)/log(2)*d+0.5)); fflush(stdout);
00433
00434 for(r=0,M=1;r<d;r++)
00435 {
00436 M=N*M;
00437 NN[r]=N;
00438 nn[r]=2*N;
00439 }
00440
00441 nfft_init_guru(&p, d, NN, M, nn, 4,
00442 PRE_PHI_HUT|
00443 PRE_PSI|
00444 MALLOC_F_HAT| MALLOC_X| MALLOC_F|
00445 FFTW_INIT| FFT_OUT_OF_PLACE,
00446 FFTW_MEASURE| FFTW_DESTROY_INPUT);
00447
00448 p_fft=fftw_plan_dft(d, NN, p.f_hat, p.f, FFTW_FORWARD, FFTW_MEASURE);
00449
00450 double _Complex *swapndft=(double _Complex*)nfft_malloc(p.M_total*sizeof(double _Complex));
00451
00453 nfft_vrand_shifted_unit_double(p.x, p.d*p.M_total);
00454
00455
00456
00457
00458
00459
00460
00461 nfft_precompute_one_psi(&p);
00462
00464 nfft_vrand_unit_complex(p.f_hat, p.N_total);
00465
00467 t_fft=0;
00468 r=0;
00469 while(t_fft<0.1)
00470 {
00471 r++;
00472 t=nfft_second();
00473 fftw_execute(p_fft);
00474 t=nfft_second()-t;
00475 t_fft+=t;
00476 }
00477 t_fft/=r;
00478
00479 printf("%.1e\t",t_fft);
00480
00482 nfft_vrand_unit_complex(p.f_hat, p.N_total);
00483
00485 if(test_ndft)
00486 {
00487 NFFT_SWAP_complex(p.f,swapndft);
00488 t_ndft=0;
00489 r=0;
00490 while(t_ndft<0.1)
00491 {
00492 r++;
00493 t=nfft_second();
00494 ndft_trafo(&p);
00495 t=nfft_second()-t;
00496 t_ndft+=t;
00497 }
00498 t_ndft/=r;
00499 printf("%.1e\t",t_ndft);
00500
00501
00502
00503 NFFT_SWAP_complex(p.f,swapndft);
00504 }
00505 else
00506 printf("\t");
00507
00509 t_nfft=0;
00510 r=0;
00511 while(t_nfft<0.1)
00512 {
00513 r++;
00514 t=nfft_second();
00515 nfft_trafo(&p);
00516 t=nfft_second()-t;
00517 t_nfft+=t;
00518 }
00519 t_nfft/=r;
00520 printf("%.1e\t",t_nfft);
00521 if(test_ndft)
00522 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f, p.M_total));
00523
00524
00525
00527 t_nfft=0;
00528 r=0;
00529 while(t_nfft<0.1)
00530 {
00531 r++;
00532 t=nfft_second();
00533 nfft_trafo_2d(&p);
00534 t=nfft_second()-t;
00535 t_nfft+=t;
00536 }
00537 t_nfft/=r;
00538 printf("%.1e\t",t_nfft);
00539 if(test_ndft)
00540 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f, p.M_total));
00541
00542
00543
00544 printf("\n");
00545
00546 nfft_free(swapndft);
00547 fftw_destroy_plan(p_fft);
00548 nfft_finalize(&p);
00549 }
00550
00551
00552 void measure_time_nfft_XXX5(int d, int N, unsigned test_ndft)
00553 {
00554 int j,r, M, NN[d], nn[d];
00555 double t, t_fft, t_ndft, t_nfft;
00556
00557 nfft_plan p;
00558 fftw_plan p_fft;
00559
00560 printf("%d\t",(int)(log(N)/log(2)*d+0.5)); fflush(stdout);
00561
00562 for(r=0,M=1;r<d;r++)
00563 {
00564 M=N*M;
00565 NN[r]=N;
00566 nn[r]=2*N;
00567 }
00568
00569 nfft_init_guru(&p, d, NN, M, nn, 4,
00570 PRE_PHI_HUT|
00571 PRE_PSI|
00572 MALLOC_F_HAT| MALLOC_X| MALLOC_F|
00573 FFTW_INIT| FFT_OUT_OF_PLACE,
00574 FFTW_MEASURE| FFTW_DESTROY_INPUT);
00575
00576 p_fft=fftw_plan_dft(d, NN, p.f, p.f_hat, FFTW_FORWARD, FFTW_MEASURE);
00577
00578 double _Complex *swapndft=(double _Complex*)nfft_malloc(p.N_total*sizeof(double _Complex));
00579
00581 nfft_vrand_shifted_unit_double(p.x, p.d*p.M_total);
00582
00583
00584
00585 nfft_precompute_one_psi(&p);
00586
00588 nfft_vrand_unit_complex(p.f, p.M_total);
00589
00591 t_fft=0;
00592 r=0;
00593 while(t_fft<0.1)
00594 {
00595 r++;
00596 t=nfft_second();
00597 fftw_execute(p_fft);
00598 t=nfft_second()-t;
00599 t_fft+=t;
00600 }
00601 t_fft/=r;
00602
00603 printf("%.1e\t",t_fft);
00604
00606 nfft_vrand_unit_complex(p.f, p.M_total);
00607
00609 if(test_ndft)
00610 {
00611 NFFT_SWAP_complex(p.f_hat,swapndft);
00612 t_ndft=0;
00613 r=0;
00614 while(t_ndft<0.1)
00615 {
00616 r++;
00617 t=nfft_second();
00618 ndft_adjoint(&p);
00619 t=nfft_second()-t;
00620 t_ndft+=t;
00621 }
00622 t_ndft/=r;
00623 printf("%.1e\t",t_ndft);
00624
00625
00626
00627 NFFT_SWAP_complex(p.f_hat,swapndft);
00628 }
00629 else
00630 printf("\t");
00631
00633 t_nfft=0;
00634 r=0;
00635 while(t_nfft<0.1)
00636 {
00637 r++;
00638 t=nfft_second();
00639 nfft_adjoint(&p);
00640 t=nfft_second()-t;
00641 t_nfft+=t;
00642 }
00643 t_nfft/=r;
00644 printf("%.1e\t",t_nfft);
00645 if(test_ndft)
00646 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f_hat, p.N_total));
00647
00648
00649
00651 t_nfft=0;
00652 r=0;
00653 while(t_nfft<0.1)
00654 {
00655 r++;
00656 t=nfft_second();
00657 nfft_adjoint_2d(&p);
00658 t=nfft_second()-t;
00659 t_nfft+=t;
00660 }
00661 t_nfft/=r;
00662 printf("%.1e\t",t_nfft);
00663 if(test_ndft)
00664 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f_hat, p.N_total));
00665
00666
00667
00668 printf("\n");
00669
00670 nfft_free(swapndft);
00671 fftw_destroy_plan(p_fft);
00672 nfft_finalize(&p);
00673 }
00674
00675
00676 void measure_time_nfft_XXX6(int d, int N, unsigned test_ndft)
00677 {
00678 int j,r, M, NN[d], nn[d];
00679 double t, t_fft, t_ndft, t_nfft;
00680
00681 nfft_plan p;
00682 fftw_plan p_fft;
00683
00684 printf("%d\t",(int)(log(N)/log(2)*d+0.5)); fflush(stdout);
00685
00686 for(r=0,M=1;r<d;r++)
00687 {
00688 M=N*M;
00689 NN[r]=N;
00690 nn[r]=2*N;
00691 }
00692
00693 nfft_init_guru(&p, d, NN, M, nn, 2,
00694 PRE_PHI_HUT|
00695 FG_PSI|
00696 MALLOC_F_HAT| MALLOC_X| MALLOC_F|
00697 FFTW_INIT| FFT_OUT_OF_PLACE,
00698 FFTW_MEASURE| FFTW_DESTROY_INPUT);
00699
00700 p_fft=fftw_plan_dft(d, NN, p.f_hat, p.f, FFTW_FORWARD, FFTW_MEASURE);
00701
00702 double _Complex *swapndft=(double _Complex*)nfft_malloc(p.M_total*sizeof(double _Complex));
00703
00705 nfft_vrand_shifted_unit_double(p.x, p.d*p.M_total);
00706
00707
00708
00709
00710 nfft_precompute_one_psi(&p);
00711
00713 nfft_vrand_unit_complex(p.f_hat, p.N_total);
00714
00716 t_fft=0;
00717 r=0;
00718 while(t_fft<0.1)
00719 {
00720 r++;
00721 t=nfft_second();
00722 fftw_execute(p_fft);
00723 t=nfft_second()-t;
00724 t_fft+=t;
00725 }
00726 t_fft/=r;
00727
00728 printf("%.1e\t",t_fft);
00729
00731 nfft_vrand_unit_complex(p.f_hat, p.N_total);
00732
00734 if(test_ndft)
00735 {
00736 NFFT_SWAP_complex(p.f,swapndft);
00737 t_ndft=0;
00738 r=0;
00739 while(t_ndft<0.1)
00740 {
00741 r++;
00742 t=nfft_second();
00743 ndft_trafo(&p);
00744 t=nfft_second()-t;
00745 t_ndft+=t;
00746 }
00747 t_ndft/=r;
00748 printf("%.1e\t",t_ndft);
00749
00750
00751
00752 NFFT_SWAP_complex(p.f,swapndft);
00753 }
00754 else
00755 printf("\t");
00756
00758 t_nfft=0;
00759 r=0;
00760 while(t_nfft<0.1)
00761 {
00762 r++;
00763 t=nfft_second();
00764 nfft_trafo(&p);
00765 t=nfft_second()-t;
00766 t_nfft+=t;
00767 }
00768 t_nfft/=r;
00769 printf("%.1e\t",t_nfft);
00770 if(test_ndft)
00771 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f, p.M_total));
00772
00773
00774
00776 t_nfft=0;
00777 r=0;
00778 while(t_nfft<0.1)
00779 {
00780 r++;
00781 t=nfft_second();
00782 nfft_trafo_3d(&p);
00783 t=nfft_second()-t;
00784 t_nfft+=t;
00785 }
00786 t_nfft/=r;
00787 printf("%.1e\t",t_nfft);
00788 if(test_ndft)
00789 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f, p.M_total));
00790
00791
00792
00793 printf("\n");
00794
00795 nfft_free(swapndft);
00796 fftw_destroy_plan(p_fft);
00797 nfft_finalize(&p);
00798 }
00799
00800
00801 void measure_time_nfft_XXX7(int d, int N, unsigned test_ndft)
00802 {
00803 int j,r, M, NN[d], nn[d];
00804 double t, t_fft, t_ndft, t_nfft;
00805
00806 nfft_plan p;
00807 fftw_plan p_fft;
00808
00809 printf("%d\t",(int)(log(N)/log(2)*d+0.5)); fflush(stdout);
00810
00811 for(r=0,M=1;r<d;r++)
00812 {
00813 M=N*M;
00814 NN[r]=N;
00815 nn[r]=2*N;
00816 }
00817
00818 nfft_init_guru(&p, d, NN, M, nn, 2,
00819 PRE_PHI_HUT|
00820 FG_PSI|
00821 MALLOC_F_HAT| MALLOC_X| MALLOC_F|
00822 FFTW_INIT| FFT_OUT_OF_PLACE,
00823 FFTW_MEASURE| FFTW_DESTROY_INPUT);
00824
00825 p_fft=fftw_plan_dft(d, NN, p.f, p.f_hat, FFTW_FORWARD, FFTW_MEASURE);
00826
00827 double _Complex *swapndft=(double _Complex*)nfft_malloc(p.N_total*sizeof(double _Complex));
00828
00830 nfft_vrand_shifted_unit_double(p.x, p.d*p.M_total);
00831
00832
00833
00834 nfft_precompute_one_psi(&p);
00835
00837 nfft_vrand_unit_complex(p.f, p.M_total);
00838
00840 t_fft=0;
00841 r=0;
00842 while(t_fft<0.1)
00843 {
00844 r++;
00845 t=nfft_second();
00846 fftw_execute(p_fft);
00847 t=nfft_second()-t;
00848 t_fft+=t;
00849 }
00850 t_fft/=r;
00851
00852 printf("%.1e\t",t_fft);
00853
00855 nfft_vrand_unit_complex(p.f, p.M_total);
00856
00858 if(test_ndft)
00859 {
00860 NFFT_SWAP_complex(p.f_hat,swapndft);
00861 t_ndft=0;
00862 r=0;
00863 while(t_ndft<0.1)
00864 {
00865 r++;
00866 t=nfft_second();
00867 ndft_adjoint(&p);
00868 t=nfft_second()-t;
00869 t_ndft+=t;
00870 }
00871 t_ndft/=r;
00872 printf("%.1e\t",t_ndft);
00873
00874
00875
00876 NFFT_SWAP_complex(p.f_hat,swapndft);
00877 }
00878 else
00879 printf("\t");
00880
00882 t_nfft=0;
00883 r=0;
00884 while(t_nfft<0.1)
00885 {
00886 r++;
00887 t=nfft_second();
00888 nfft_adjoint(&p);
00889 t=nfft_second()-t;
00890 t_nfft+=t;
00891 }
00892 t_nfft/=r;
00893 printf("%.1e\t",t_nfft);
00894 if(test_ndft)
00895 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f_hat, p.N_total));
00896
00897
00898
00900 t_nfft=0;
00901 r=0;
00902 while(t_nfft<0.1)
00903 {
00904 r++;
00905 t=nfft_second();
00906 nfft_adjoint_3d(&p);
00907 t=nfft_second()-t;
00908 t_nfft+=t;
00909 }
00910 t_nfft/=r;
00911 printf("%.1e\t",t_nfft);
00912 if(test_ndft)
00913 printf("(%.1e)\t",nfft_error_l_2_complex(swapndft, p.f_hat, p.N_total));
00914
00915
00916
00917 printf("\n");
00918
00919 nfft_free(swapndft);
00920 fftw_destroy_plan(p_fft);
00921 nfft_finalize(&p);
00922 }
00923
00924 int main2(void)
00925 {
00926 int l,d,logIN;
00927
00928 for(l=3;l<=6;l++)
00929 {
00930 d=3;
00931 logIN=d*l;
00932 if(logIN<=15)
00933 {
00934 measure_time_nfft_XXX6(d,(1U<< (logIN/d)),1);
00935 measure_time_nfft_XXX7(d,(1U<< (logIN/d)),1);
00936 }
00937 else
00938 {
00939 measure_time_nfft_XXX6(d,(1U<< (logIN/d)),0);
00940 measure_time_nfft_XXX7(d,(1U<< (logIN/d)),0);
00941 }
00942 }
00943
00944 exit(-1);
00945
00946
00947 for(l=7;l<=12;l++)
00948 {
00949 d=2;
00950 logIN=d*l;
00951 if(logIN<=15)
00952 {
00953 measure_time_nfft_XXX4(d,(1U<< (logIN/d)),1);
00954 measure_time_nfft_XXX5(d,(1U<< (logIN/d)),1);
00955 }
00956 else
00957 {
00958 measure_time_nfft_XXX4(d,(1U<< (logIN/d)),0);
00959 measure_time_nfft_XXX5(d,(1U<< (logIN/d)),0);
00960 }
00961 }
00962
00963 exit(-1);
00964
00965 for(l=3;l<=12;l++)
00966 {
00967 logIN=l;
00968 if(logIN<=15)
00969 {
00970 measure_time_nfft_XXX2(1,(1U<< (logIN)),1);
00971 measure_time_nfft_XXX3(1,(1U<< (logIN)),1);
00972 }
00973 else
00974 {
00975 measure_time_nfft_XXX2(1,(1U<< (logIN)),0);
00976 measure_time_nfft_XXX3(1,(1U<< (logIN)),0);
00977 }
00978 }
00979
00980 exit(-1);
00981 }
00982
00983 int XXX(void)
00984 {
00985 int l,d,logIN;
00986
00987 }
00988
00989 int main(void)
00990 {
00991 int l,d,logIN;
00992
00993 printf("\\hline $l_N$ & FFT & NDFT & NFFT & NFFT/FFT\\\\\n");
00994 printf("\\hline \\hline \\multicolumn{5}{|c|}{$d=1$} \\\\ \\hline\n");
00995 for(l=3;l<=22;l++)
00996 {
00997 d=1;
00998 logIN=l;
00999 if(logIN<=15)
01000 measure_time_nfft(d,(1U<< (logIN/d)),1);
01001 else
01002 measure_time_nfft(d,(1U<< (logIN/d)),0);
01003
01004 fflush(stdout);
01005 }
01006
01007 printf("\\hline $l_N$ & FFT & NDFT & NFFT & NFFT/FFT\\\\\n");
01008 printf("\\hline \\hline \\multicolumn{5}{|c|}{$d=2$} \\\\ \\hline\n");
01009 for(l=3;l<=11;l++)
01010 {
01011 d=2;
01012 logIN=d*l;
01013 if(logIN<=15)
01014 measure_time_nfft(d,(1U<< (logIN/d)),1);
01015 else
01016 measure_time_nfft(d,(1U<< (logIN/d)),0);
01017
01018 fflush(stdout);
01019 }
01020
01021 printf("\\hline \\hline \\multicolumn{5}{|c|}{$d=3$} \\\\ \\hline\n");
01022 for(l=3;l<=7;l++)
01023 {
01024 d=3;
01025 logIN=d*l;
01026 if(logIN<=15)
01027 measure_time_nfft(d,(1U<< (logIN/d)),1);
01028 else
01029 measure_time_nfft(d,(1U<< (logIN/d)),0);
01030
01031 fflush(stdout);
01032 }
01033
01034 return 1;
01035 }