46 enum boolean {NO = 0, YES = 1};
52 enum gridtype {GRID_GAUSS_LEGENDRE = 0, GRID_CLENSHAW_CURTIS = 1,
53 GRID_HEALPIX = 2, GRID_EQUIDISTRIBUTION = 3, GRID_EQUIDISTRIBUTION_UNIFORM = 4};
56 enum functiontype {FUNCTION_RANDOM_BANDLIMITED = 0, FUNCTION_F1 = 1,
57 FUNCTION_F2 = 2, FUNCTION_F3 = 3, FUNCTION_F4 = 4, FUNCTION_F5 = 5,
61 enum modes {USE_GRID = 0, RANDOM = 1};
71 int main (
int argc,
char **argv)
103 double _Complex *f_grid;
104 double _Complex *f_compare;
106 double _Complex *f_hat_gen;
108 double _Complex *f_hat;
115 double err_infty_avg;
138 double x1,x2,x3,temp;
147 fscanf(stdin,
"testcases=%d\n",&tc_max);
148 fprintf(stdout,
"%d\n",tc_max);
151 for (tc = 0; tc < tc_max; tc++)
154 fscanf(stdin,
"nfsft=%d\n",&use_nfsft);
155 fprintf(stdout,
"%d\n",use_nfsft);
159 fscanf(stdin,
"nfft=%d\n",&use_nfft);
160 fprintf(stdout,
"%d\n",use_nfsft);
164 fscanf(stdin,
"cutoff=%d\n",&cutoff);
165 fprintf(stdout,
"%d\n",cutoff);
174 fscanf(stdin,
"fpt=%d\n",&use_fpt);
175 fprintf(stdout,
"%d\n",use_fpt);
179 fscanf(stdin,
"threshold=%lf\n",&threshold);
180 fprintf(stdout,
"%lf\n",threshold);
200 fscanf(stdin,
"testmode=%d\n",&testmode);
201 fprintf(stdout,
"%d\n",testmode);
203 if (testmode == ERROR)
206 fscanf(stdin,
"gridtype=%d\n",&gridtype);
207 fprintf(stdout,
"%d\n",gridtype);
210 fscanf(stdin,
"testfunction=%d\n",&testfunction);
211 fprintf(stdout,
"%d\n",testfunction);
214 if (testfunction == FUNCTION_RANDOM_BANDLIMITED)
217 fscanf(stdin,
"bandlimit=%d\n",&N);
218 fprintf(stdout,
"%d\n",N);
226 fscanf(stdin,
"repetitions=%d\n",&repetitions);
227 fprintf(stdout,
"%d\n",repetitions);
229 fscanf(stdin,
"mode=%d\n",&mode);
230 fprintf(stdout,
"%d\n",mode);
235 fscanf(stdin,
"points=%d\n",&m_compare);
236 fprintf(stdout,
"%d\n",m_compare);
237 x_compare = (
double*)
nfft_malloc(2*m_compare*
sizeof(
double));
239 while (d < m_compare)
241 x1 = 2.0*(((double)rand())/RAND_MAX) - 1.0;
242 x2 = 2.0*(((double)rand())/RAND_MAX) - 1.0;
243 x3 = 2.0*(((double)rand())/RAND_MAX) - 1.0;
244 temp = sqrt(x1*x1+x2*x2+x3*x3);
247 x_compare[2*d+1] = acos(x3);
248 if (x_compare[2*d+1] == 0 || x_compare[2*d+1] == KPI)
250 x_compare[2*d] = 0.0;
254 x_compare[2*d] = atan2(x2/sin(x_compare[2*d+1]),x1/sin(x_compare[2*d+1]));
256 x_compare[2*d] *= 1.0/(2.0*KPI);
257 x_compare[2*d+1] *= 1.0/(2.0*KPI);
261 f_compare = (
double _Complex*)
nfft_malloc(m_compare*
sizeof(
double _Complex));
262 f = (
double _Complex*)
nfft_malloc(m_compare*
sizeof(
double _Complex));
271 fscanf(stdin,
"bandwidths=%d\n",&iNQ_max);
272 fprintf(stdout,
"%d\n",iNQ_max);
277 if (testmode == TIMING)
283 for (iNQ = 0; iNQ < iNQ_max; iNQ++)
285 if (testmode == TIMING)
288 fscanf(stdin,
"%d %d %d\n",&NQ[iNQ],&SQ[iNQ],&RQ[iNQ]);
289 fprintf(stdout,
"%d %d %d\n",NQ[iNQ],SQ[iNQ],RQ[iNQ]);
290 NQ_max = MAX(NQ_max,NQ[iNQ]);
291 SQ_max = MAX(SQ_max,SQ[iNQ]);
296 fscanf(stdin,
"%d %d\n",&NQ[iNQ],&SQ[iNQ]);
297 fprintf(stdout,
"%d %d\n",NQ[iNQ],SQ[iNQ]);
298 NQ_max = MAX(NQ_max,NQ[iNQ]);
299 SQ_max = MAX(SQ_max,SQ[iNQ]);
306 nfsft_precompute(NQ_max, threshold,
307 ((use_nfsft==NO)?(NFSFT_NO_FAST_ALGORITHM):(0U)), 0U);
309 if (testmode == TIMING)
312 f_hat = (
double _Complex*)
nfft_malloc(NFSFT_F_HAT_SIZE(NQ_max)*
sizeof(
double _Complex));
313 f = (
double _Complex*)
nfft_malloc(SQ_max*
sizeof(
double _Complex));
314 x_grid = (
double*)
nfft_malloc(2*SQ_max*
sizeof(
double));
315 for (d = 0; d < SQ_max; d++)
317 f[d] = (((double)rand())/RAND_MAX)-0.5 + _Complex_I*((((double)rand())/RAND_MAX)-0.5);
318 x_grid[2*d] = (((double)rand())/RAND_MAX) - 0.5;
319 x_grid[2*d+1] = (((double)rand())/RAND_MAX) * 0.5;
326 for (iNQ = 0; iNQ < iNQ_max; iNQ++)
328 if (testmode == TIMING)
330 nfsft_init_guru(&plan,NQ[iNQ],SQ[iNQ], NFSFT_NORMALIZED |
331 ((use_nfft!=NO)?(0U):(NFSFT_USE_NDFT)) |
332 ((use_fpt!=NO)?(0U):(NFSFT_USE_DPT)),
333 PRE_PHI_HUT | PRE_PSI | FFTW_INIT | FFTW_MEASURE | FFT_OUT_OF_PLACE,
340 nfsft_precompute_x(&plan);
344 for (i = 0; i < RQ[iNQ]; i++)
351 nfsft_adjoint(&plan);
356 nfsft_adjoint_direct(&plan);
363 t_avg = t_avg/((double)RQ[iNQ]);
365 nfsft_finalize(&plan);
367 fprintf(stdout,
"%+le\n", t_avg);
368 fprintf(stderr,
"%d: %4d %4d %+le\n", tc, NQ[iNQ], SQ[iNQ], t_avg);
375 case GRID_GAUSS_LEGENDRE:
377 m_theta = SQ[iNQ] + 1;
378 m_phi = 2*SQ[iNQ] + 2;
379 m_total = m_theta*m_phi;
381 case GRID_CLENSHAW_CURTIS:
383 m_theta = 2*SQ[iNQ] + 1;
384 m_phi = 2*SQ[iNQ] + 2;
385 m_total = m_theta*m_phi;
389 m_phi = 12*SQ[iNQ]*SQ[iNQ];
390 m_total = m_theta * m_phi;
393 case GRID_EQUIDISTRIBUTION:
394 case GRID_EQUIDISTRIBUTION_UNIFORM:
397 for (k = 1; k < SQ[iNQ]; k++)
399 m_theta += (int)floor((2*KPI)/acos((cos(KPI/(
double)SQ[iNQ])-
400 cos(k*KPI/(
double)SQ[iNQ])*cos(k*KPI/(
double)SQ[iNQ]))/
401 (sin(k*KPI/(
double)SQ[iNQ])*sin(k*KPI/(
double)SQ[iNQ]))));
406 m_total = m_theta * m_phi;
412 x_grid = (
double*)
nfft_malloc(2*m_total*
sizeof(
double));
418 case GRID_GAUSS_LEGENDRE:
423 for (k = 0; k < m_theta; k++)
425 fscanf(stdin,
"%le\n",&w[k]);
426 w[k] *= (2.0*KPI)/((
double)m_phi);
432 theta = (
double*)
nfft_malloc(m_theta*
sizeof(
double));
443 for (k = 0; k < m_theta; k++)
445 fscanf(stdin,
"%le\n",&theta[k]);
449 for (n = 0; n < m_phi; n++)
451 phi[n] = n/((double)m_phi);
452 phi[n] -= ((phi[n]>=0.5)?(1.0):(0.0));
460 for (k = 0; k < m_theta; k++)
462 for (n = 0; n < m_phi; n++)
464 x_grid[2*d] = phi[n];
465 x_grid[2*d+1] = theta[k];
478 case GRID_CLENSHAW_CURTIS:
481 theta = (
double*)
nfft_malloc(m_theta*
sizeof(
double));
485 for (k = 0; k < m_theta; k++)
487 theta[k] = k/((double)2*(m_theta-1));
491 for (n = 0; n < m_phi; n++)
493 phi[n] = n/((double)m_phi);
494 phi[n] -= ((phi[n]>=0.5)?(1.0):(0.0));
498 fplan = fftw_plan_r2r_1d(SQ[iNQ]+1, w, w, FFTW_REDFT00, 0U);
499 for (k = 0; k < SQ[iNQ]+1; k++)
501 w[k] = -2.0/(4*k*k-1);
506 for (k = 0; k < SQ[iNQ]+1; k++)
508 w[k] *= (2.0*KPI)/((
double)(m_theta-1)*m_phi);
509 w[m_theta-1-k] = w[k];
511 fftw_destroy_plan(fplan);
515 for (k = 0; k < m_theta; k++)
517 for (n = 0; n < m_phi; n++)
519 x_grid[2*d] = phi[n];
520 x_grid[2*d+1] = theta[k];
534 for (k = 1; k <= SQ[iNQ]-1; k++)
536 for (n = 0; n <= 4*k-1; n++)
538 x_grid[2*d+1] = 1 - (k*k)/((
double)(3.0*SQ[iNQ]*SQ[iNQ]));
539 x_grid[2*d] = ((n+0.5)/(4*k));
540 x_grid[2*d] -= (x_grid[2*d]>=0.5)?(1.0):(0.0);
547 for (k = SQ[iNQ]; k <= 3*SQ[iNQ]; k++)
549 for (n = 0; n <= 4*SQ[iNQ]-1; n++)
551 x_grid[2*d+1] = 2.0/(3*SQ[iNQ])*(2*SQ[iNQ]-k);
552 x_grid[2*d] = (n+((k%2==0)?(0.5):(0.0)))/(4*SQ[iNQ]);
553 x_grid[2*d] -= (x_grid[2*d]>=0.5)?(1.0):(0.0);
558 for (k = 1; k <= SQ[iNQ]-1; k++)
560 for (n = 0; n <= 4*k-1; n++)
562 x_grid[2*d+1] = -x_grid[2*d2+1];
563 x_grid[2*d] = x_grid[2*d2];
569 for (d = 0; d < m_total; d++)
571 x_grid[2*d+1] = acos(x_grid[2*d+1])/(2.0*KPI);
574 w[0] = (4.0*KPI)/(m_total);
577 case GRID_EQUIDISTRIBUTION:
578 case GRID_EQUIDISTRIBUTION_UNIFORM:
581 if (gridtype == GRID_EQUIDISTRIBUTION)
583 w_temp = (
double*)
nfft_malloc((SQ[iNQ]+1)*
sizeof(double));
584 fplan = fftw_plan_r2r_1d(SQ[iNQ]/2+1, w_temp, w_temp, FFTW_REDFT00, 0U);
585 for (k = 0; k < SQ[iNQ]/2+1; k++)
587 w_temp[k] = -2.0/(4*k*k-1);
592 for (k = 0; k < SQ[iNQ]/2+1; k++)
594 w_temp[k] *= (2.0*KPI)/((
double)(SQ[iNQ]));
595 w_temp[SQ[iNQ]-k] = w_temp[k];
597 fftw_destroy_plan(fplan);
603 if (gridtype == GRID_EQUIDISTRIBUTION)
609 w[d] = (4.0*KPI)/(m_total);
614 if (gridtype == GRID_EQUIDISTRIBUTION)
616 w[d] = w_temp[SQ[iNQ]];
620 w[d] = (4.0*KPI)/(m_total);
624 for (k = 1; k < SQ[iNQ]; k++)
626 theta_s = (double)k*KPI/(
double)SQ[iNQ];
627 M = (int)floor((2.0*KPI)/acos((cos(KPI/(
double)SQ[iNQ])-
628 cos(theta_s)*cos(theta_s))/(sin(theta_s)*sin(theta_s))));
630 for (n = 0; n < M; n++)
632 x_grid[2*d] = (n + 0.5)/M;
633 x_grid[2*d] -= (x_grid[2*d]>=0.5)?(1.0):(0.0);
634 x_grid[2*d+1] = theta_s/(2.0*KPI);
635 if (gridtype == GRID_EQUIDISTRIBUTION)
637 w[d] = w_temp[k]/((double)(M));
641 w[d] = (4.0*KPI)/(m_total);
647 if (gridtype == GRID_EQUIDISTRIBUTION)
658 f_grid = (
double _Complex*)
nfft_malloc(m_total*
sizeof(
double _Complex));
666 f_compare = (
double _Complex*)
nfft_malloc(m_compare*
sizeof(
double _Complex));
673 switch (testfunction)
675 case FUNCTION_RANDOM_BANDLIMITED:
676 f_hat_gen = (
double _Complex*)
nfft_malloc(NFSFT_F_HAT_SIZE(N)*
sizeof(
double _Complex));
681 nfsft_init_guru(&plan_gen,N,m_total, NFSFT_NORMALIZED |
682 ((use_nfft!=NO)?(0U):(NFSFT_USE_NDFT)) |
683 ((use_fpt!=NO)?(0U):(NFSFT_USE_DPT)),
684 ((N>512)?(0U):(PRE_PHI_HUT | PRE_PSI)) | FFTW_INIT |
685 FFT_OUT_OF_PLACE, cutoff);
687 plan_gen.
f_hat = f_hat_gen;
691 nfsft_precompute_x(&plan_gen);
693 for (k = 0; k < plan_gen.
N_total; k++)
698 for (k = 0; k <= N; k++)
700 for (n = -k; n <= k; n++)
702 f_hat_gen[NFSFT_INDEX(k,n,&plan_gen)] =
703 (((double)rand())/RAND_MAX)-0.5 + _Complex_I*((((double)rand())/RAND_MAX)-0.5);
710 nfsft_trafo(&plan_gen);
715 nfsft_trafo_direct(&plan_gen);
718 nfsft_finalize(&plan_gen);
722 nfsft_init_guru(&plan_gen,N,m_compare, NFSFT_NORMALIZED |
723 ((use_nfft!=NO)?(0U):(NFSFT_USE_NDFT)) |
724 ((use_fpt!=NO)?(0U):(NFSFT_USE_DPT)),
725 ((N>512)?(0U):(PRE_PHI_HUT | PRE_PSI)) | FFTW_INIT |
726 FFT_OUT_OF_PLACE, cutoff);
728 plan_gen.
f_hat = f_hat_gen;
729 plan_gen.
x = x_compare;
730 plan_gen.
f = f_compare;
732 nfsft_precompute_x(&plan_gen);
737 nfsft_trafo(&plan_gen);
742 nfsft_trafo_direct(&plan_gen);
745 nfsft_finalize(&plan_gen);
749 memcpy(f_compare,f_grid,m_total*
sizeof(
double _Complex));
757 for (d = 0; d < m_total; d++)
759 x1 = sin(x_grid[2*d+1]*2.0*KPI)*cos(x_grid[2*d]*2.0*KPI);
760 x2 = sin(x_grid[2*d+1]*2.0*KPI)*sin(x_grid[2*d]*2.0*KPI);
761 x3 = cos(x_grid[2*d+1]*2.0*KPI);
762 f_grid[d] = x1*x2*x3;
766 for (d = 0; d < m_compare; d++)
768 x1 = sin(x_compare[2*d+1]*2.0*KPI)*cos(x_compare[2*d]*2.0*KPI);
769 x2 = sin(x_compare[2*d+1]*2.0*KPI)*sin(x_compare[2*d]*2.0*KPI);
770 x3 = cos(x_compare[2*d+1]*2.0*KPI);
771 f_compare[d] = x1*x2*x3;
776 memcpy(f_compare,f_grid,m_total*
sizeof(
double _Complex));
780 for (d = 0; d < m_total; d++)
782 x1 = sin(x_grid[2*d+1]*2.0*KPI)*cos(x_grid[2*d]*2.0*KPI);
783 x2 = sin(x_grid[2*d+1]*2.0*KPI)*sin(x_grid[2*d]*2.0*KPI);
784 x3 = cos(x_grid[2*d+1]*2.0*KPI);
785 f_grid[d] = 0.1*exp(x1+x2+x3);
789 for (d = 0; d < m_compare; d++)
791 x1 = sin(x_compare[2*d+1]*2.0*KPI)*cos(x_compare[2*d]*2.0*KPI);
792 x2 = sin(x_compare[2*d+1]*2.0*KPI)*sin(x_compare[2*d]*2.0*KPI);
793 x3 = cos(x_compare[2*d+1]*2.0*KPI);
794 f_compare[d] = 0.1*exp(x1+x2+x3);
799 memcpy(f_compare,f_grid,m_total*
sizeof(
double _Complex));
803 for (d = 0; d < m_total; d++)
805 x1 = sin(x_grid[2*d+1]*2.0*KPI)*cos(x_grid[2*d]*2.0*KPI);
806 x2 = sin(x_grid[2*d+1]*2.0*KPI)*sin(x_grid[2*d]*2.0*KPI);
807 x3 = cos(x_grid[2*d+1]*2.0*KPI);
808 temp = sqrt(x1*x1)+sqrt(x2*x2)+sqrt(x3*x3);
809 f_grid[d] = 0.1*temp;
813 for (d = 0; d < m_compare; d++)
815 x1 = sin(x_compare[2*d+1]*2.0*KPI)*cos(x_compare[2*d]*2.0*KPI);
816 x2 = sin(x_compare[2*d+1]*2.0*KPI)*sin(x_compare[2*d]*2.0*KPI);
817 x3 = cos(x_compare[2*d+1]*2.0*KPI);
818 temp = sqrt(x1*x1)+sqrt(x2*x2)+sqrt(x3*x3);
819 f_compare[d] = 0.1*temp;
824 memcpy(f_compare,f_grid,m_total*
sizeof(
double _Complex));
828 for (d = 0; d < m_total; d++)
830 x1 = sin(x_grid[2*d+1]*2.0*KPI)*cos(x_grid[2*d]*2.0*KPI);
831 x2 = sin(x_grid[2*d+1]*2.0*KPI)*sin(x_grid[2*d]*2.0*KPI);
832 x3 = cos(x_grid[2*d+1]*2.0*KPI);
833 temp = sqrt(x1*x1)+sqrt(x2*x2)+sqrt(x3*x3);
834 f_grid[d] = 1.0/(temp);
838 for (d = 0; d < m_compare; d++)
840 x1 = sin(x_compare[2*d+1]*2.0*KPI)*cos(x_compare[2*d]*2.0*KPI);
841 x2 = sin(x_compare[2*d+1]*2.0*KPI)*sin(x_compare[2*d]*2.0*KPI);
842 x3 = cos(x_compare[2*d+1]*2.0*KPI);
843 temp = sqrt(x1*x1)+sqrt(x2*x2)+sqrt(x3*x3);
844 f_compare[d] = 1.0/(temp);
849 memcpy(f_compare,f_grid,m_total*
sizeof(
double _Complex));
853 for (d = 0; d < m_total; d++)
855 x1 = sin(x_grid[2*d+1]*2.0*KPI)*cos(x_grid[2*d]*2.0*KPI);
856 x2 = sin(x_grid[2*d+1]*2.0*KPI)*sin(x_grid[2*d]*2.0*KPI);
857 x3 = cos(x_grid[2*d+1]*2.0*KPI);
858 temp = sqrt(x1*x1)+sqrt(x2*x2)+sqrt(x3*x3);
859 f_grid[d] = 0.1*sin(1+temp)*sin(1+temp);
863 for (d = 0; d < m_compare; d++)
865 x1 = sin(x_compare[2*d+1]*2.0*KPI)*cos(x_compare[2*d]*2.0*KPI);
866 x2 = sin(x_compare[2*d+1]*2.0*KPI)*sin(x_compare[2*d]*2.0*KPI);
867 x3 = cos(x_compare[2*d+1]*2.0*KPI);
868 temp = sqrt(x1*x1)+sqrt(x2*x2)+sqrt(x3*x3);
869 f_compare[d] = 0.1*sin(1+temp)*sin(1+temp);
874 memcpy(f_compare,f_grid,m_total*
sizeof(
double _Complex));
878 for (d = 0; d < m_total; d++)
880 if (x_grid[2*d+1] <= 0.25)
886 f_grid[d] = 1.0/(sqrt(1+3*cos(2.0*KPI*x_grid[2*d+1])*cos(2.0*KPI*x_grid[2*d+1])));
891 for (d = 0; d < m_compare; d++)
893 if (x_compare[2*d+1] <= 0.25)
899 f_compare[d] = 1.0/(sqrt(1+3*cos(2.0*KPI*x_compare[2*d+1])*cos(2.0*KPI*x_compare[2*d+1])));
905 memcpy(f_compare,f_grid,m_total*
sizeof(
double _Complex));
911 for (d = 0; d < m_total; d++)
917 for (d = 0; d < m_compare; d++)
924 memcpy(f_compare,f_grid,m_total*
sizeof(
double _Complex));
932 nfsft_init_guru(&plan_adjoint,NQ[iNQ],m_total, NFSFT_NORMALIZED |
933 ((use_nfft!=NO)?(0U):(NFSFT_USE_NDFT)) |
934 ((use_fpt!=NO)?(0U):(NFSFT_USE_DPT)),
935 ((NQ[iNQ]>512)?(0U):(PRE_PHI_HUT | PRE_PSI)) | FFTW_INIT |
936 FFT_OUT_OF_PLACE, cutoff);
938 plan_adjoint_ptr = &plan_adjoint;
942 nfsft_init_guru(&plan,NQ[iNQ],m_compare, NFSFT_NORMALIZED |
943 ((use_nfft!=NO)?(0U):(NFSFT_USE_NDFT)) |
944 ((use_fpt!=NO)?(0U):(NFSFT_USE_DPT)),
945 ((NQ[iNQ]>512)?(0U):(PRE_PHI_HUT | PRE_PSI)) | FFTW_INIT |
946 FFT_OUT_OF_PLACE, cutoff);
951 plan_ptr = &plan_adjoint;
954 f_hat = (
double _Complex*)
nfft_malloc(NFSFT_F_HAT_SIZE(NQ[iNQ])*
sizeof(
double _Complex));
956 plan_adjoint_ptr->
f_hat = f_hat;
957 plan_adjoint_ptr->
x = x_grid;
958 plan_adjoint_ptr->
f = f_grid;
960 plan_ptr->
f_hat = f_hat;
961 plan_ptr->
x = x_compare;
966 nfsft_precompute_x(plan_adjoint_ptr);
967 if (plan_adjoint_ptr != plan_ptr)
969 nfsft_precompute_x(plan_ptr);
978 for (i = 0; i < 1; i++)
993 for (k = 0; k < m_theta; k++)
995 for (n = 0; n < m_phi; n++)
1023 if (use_nfsft != NO)
1026 nfsft_adjoint(plan_adjoint_ptr);
1031 nfsft_adjoint_direct(plan_adjoint_ptr);
1047 if (use_nfsft != NO)
1050 nfsft_trafo(plan_ptr);
1055 nfsft_trafo_direct(plan_ptr);
1064 nfsft_finalize(plan_adjoint_ptr);
1065 if (plan_ptr != plan_adjoint_ptr)
1067 nfsft_finalize(plan_ptr);
1074 err_infty_avg +=
X(error_l_infty_complex)(f, f_compare, m_compare);
1075 err_2_avg +=
X(error_l_2_complex)(f, f_compare, m_compare);
1097 t_avg = t_avg/((double)repetitions);
1100 err_infty_avg = err_infty_avg/((double)repetitions);
1103 err_2_avg = err_2_avg/((double)repetitions);
1106 fprintf(stdout,
"%+le %+le %+le\n", t_avg, err_infty_avg, err_2_avg);
1107 fprintf(stderr,
"%d: %4d %4d %+le %+le %+le\n", tc, NQ[iNQ], SQ[iNQ],
1108 t_avg, err_infty_avg, err_2_avg);
1112 fprintf(stderr,
"\n");
1120 if (testmode == TIMING)
1132 if (testmode == TIMING)
1143 return EXIT_SUCCESS;
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...
modes
TODO Add comment here.
testtype
Enumeration for test modes.
#define X(name)
Include header for C99 complex datatype.
functiontype
Enumeration for test functions.
void * nfft_malloc(size_t n)
NFFT_INT N_total
Total number of Fourier coefficients.
double threshold
The threshold /f$/f$.
gridtype
Enumeration for quadrature grid types.