46 fpt_set set =
fpt_init(1,lrint(ceil(log2((
double)N))),0U);
49 double *alpha = malloc((N+2)*
sizeof(
double)),
50 *beta = malloc((N+2)*
sizeof(
double)),
51 *gamma = malloc((N+2)*
sizeof(
double));
54 alpha[0] = beta[0] = 0.0;
61 for (k = 0; k <= N; k++)
63 alpha[k+1] = ((double)(2*k+1))/((
double)(k+1));
65 gamma[k+1] = -((double)(k))/((
double)(k+1));
70 "Computing a fast polynomial transform (FPT) and a fast discrete cosine \n"
71 "transform (DCT) to evaluate\n\n"
72 " f_j = a_0 P_0(x_j) + a_1 P_1(x_j) + ... + a_N P_N(x_j), j=0,1,...,N,\n\n"
73 "with N=%d, x_j = cos(j*pi/N), j=0,1,...N, the Chebyshev nodes, a_k,\n"
74 "k=0,1,...,N, random Fourier coefficients in [-1,1]x[-1,1]*I, and P_k,\n"
75 "k=0,1,...,N, the Legendre polynomials.",N
93 double _Complex *a = malloc((N+1)*
sizeof(
double _Complex));
94 double _Complex *b = malloc((N+1)*
sizeof(
double _Complex));
95 double *f = malloc((N+1)*
sizeof(
double _Complex));
98 const int NP1 = N + 1;
99 fftw_r2r_kind kind = FFTW_REDFT00;
100 fftw_plan p = fftw_plan_many_r2r(1, &NP1, 1, (
double*)b, NULL, 2, 1,
101 (
double*)f, NULL, 1, 1, &kind, 0U);
106 printf(
"\n2) Random Fourier coefficients a_k, k=0,1,...,N:\n");
107 for (k = 0; k <= N; k++)
109 a[k] = 2.0*nfft_drand48() - 1.0;
110 printf(
" a_%-2d = %+5.3lE\n",k,creal(a[k]));
115 fpt_trafo(set,0,a,b,N,0U);
123 for (j = 1; j < N; j++)
132 printf(
"\n3) Function values f_j, j=1,1,...,M:\n");
133 for (j = 0; j <= N; j++)
134 printf(
" f_%-2d = %+5.3lE\n",j,f[j]);
143 fftw_destroy_plan(p);
fpt_set fpt_init(const int M, const int t, const unsigned int flags)
void fpt_precompute(fpt_set set, const int m, double *alpha, double *beta, double *gam, int k_start, const double threshold)
Header file for the nfft3 library.