jianmu
/
lammps
forked from lijiext/lammps
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge pull request #372 from akohlmey/fft-cleanup 

simplify FFT3d code by removing support for outdated FFT libraries
This commit is contained in:
sjplimp 2017-02-13 08:53:13 -07:00 committed by GitHub
parent 2bdda8f6c0 191453e1c7
commit 4843296d4e
18 changed files with 39 additions and 608 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
doc/src/Section_start.txt
View File

@ -413,7 +413,7 @@ uses (for performing 1d FFTs) when running the particle-particle
particle-mesh (PPPM) option for long-range Coulombics via the
"kspace_style"_kspace_style.html command.
LAMMPS supports various open-source or vendor-supplied FFT libraries
LAMMPS supports common open-source or vendor-supplied FFT libraries
for this purpose. If you leave these 3 variables blank, LAMMPS will
use the open-source "KISS FFT library"_http://kissfft.sf.net, which is
included in the LAMMPS distribution. This library is portable to all
@ -423,10 +423,9 @@ package in your build, you can also leave the 3 variables blank.
Otherwise, select which kinds of FFTs to use as part of the FFT_INC
setting by a switch of the form -DFFT_XXX. Recommended values for XXX
are: MKL, SCSL, FFTW2, and FFTW3. Legacy options are: INTEL, SGI,
ACML, and T3E. For backward compatability, using -DFFT_FFTW will use
the FFTW2 library. Using -DFFT_NONE will use the KISS library
described above.
are: MKL or FFTW3. FFTW2 and NONE are supported as legacy options.
Selecting -DFFT_FFTW will use the FFTW3 library and -DFFT_NONE will
use the KISS library described above.
You may also need to set the FFT_INC, FFT_PATH, and FFT_LIB variables,
so the compiler and linker can find the needed FFT header and library

450
src/KSPACE/fft3d.cpp
View File

@ -75,20 +75,7 @@ void fft_3d(FFT_DATA *in, FFT_DATA *out, int flag, struct fft_plan_3d *plan)
// system specific constants
#if defined(FFT_SCSL)
int isys = 0;
FFT_PREC scalef = 1.0;
#elif defined(FFT_DEC)
char c = 'C';
char f = 'F';
char b = 'B';
int one = 1;
#elif defined(FFT_T3E)
int isys = 0;
double scalef = 1.0;
#elif defined(FFT_ACML)
int info;
#elif defined(FFT_FFTW3)
#if defined(FFT_FFTW3)
FFTW_API(plan) theplan;
#else
// nothing to do for other FFTs.
@ -112,35 +99,11 @@ void fft_3d(FFT_DATA *in, FFT_DATA *out, int flag, struct fft_plan_3d *plan)
total = plan->total1;
length = plan->length1;
#if defined(FFT_SGI)
for (offset = 0; offset < total; offset += length)
FFT_1D(flag,length,&data[offset],1,plan->coeff1);
#elif defined(FFT_SCSL)
for (offset = 0; offset < total; offset += length)
FFT_1D(flag,length,scalef,&data[offset],&data[offset],plan->coeff1,
plan->work1,&isys);
#elif defined(FFT_ACML)
num=total/length;
FFT_1D(&flag,&num,&length,data,plan->coeff1,&info);
#elif defined(FFT_INTEL)
for (offset = 0; offset < total; offset += length)
FFT_1D(&data[offset],&length,&flag,plan->coeff1);
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
if (flag == -1)
DftiComputeForward(plan->handle_fast,data);
else
DftiComputeBackward(plan->handle_fast,data);
#elif defined(FFT_DEC)
if (flag == -1)
for (offset = 0; offset < total; offset += length)
FFT_1D(&c,&c,&f,&data[offset],&data[offset],&length,&one);
else
for (offset = 0; offset < total; offset += length)
FFT_1D(&c,&c,&b,&data[offset],&data[offset],&length,&one);
#elif defined(FFT_T3E)
for (offset = 0; offset < total; offset += length)
FFT_1D(&flag,&length,&scalef,&data[offset],&data[offset],plan->coeff1,
plan->work1,&isys);
#elif defined(FFT_FFTW2)
if (flag == -1)
fftw(plan->plan_fast_forward,total/length,data,1,length,NULL,0,0);
@ -175,35 +138,11 @@ void fft_3d(FFT_DATA *in, FFT_DATA *out, int flag, struct fft_plan_3d *plan)
total = plan->total2;
length = plan->length2;
#if defined(FFT_SGI)
for (offset = 0; offset < total; offset += length)
FFT_1D(flag,length,&data[offset],1,plan->coeff2);
#elif defined(FFT_SCSL)
for (offset = 0; offset < total; offset += length)
FFT_1D(flag,length,scalef,&data[offset],&data[offset],plan->coeff2,
plan->work2,&isys);
#elif defined(FFT_ACML)
num=total/length;
FFT_1D(&flag,&num,&length,data,plan->coeff2,&info);
#elif defined(FFT_INTEL)
for (offset = 0; offset < total; offset += length)
FFT_1D(&data[offset],&length,&flag,plan->coeff2);
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
if (flag == -1)
DftiComputeForward(plan->handle_mid,data);
else
DftiComputeBackward(plan->handle_mid,data);
#elif defined(FFT_DEC)
if (flag == -1)
for (offset = 0; offset < total; offset += length)
FFT_1D(&c,&c,&f,&data[offset],&data[offset],&length,&one);
else
for (offset = 0; offset < total; offset += length)
FFT_1D(&c,&c,&b,&data[offset],&data[offset],&length,&one);
#elif defined(FFT_T3E)
for (offset = 0; offset < total; offset += length)
FFT_1D(&flag,&length,&scalef,&data[offset],&data[offset],plan->coeff2,
plan->work2,&isys);
#elif defined(FFT_FFTW2)
if (flag == -1)
fftw(plan->plan_mid_forward,total/length,data,1,length,NULL,0,0);
@ -238,35 +177,11 @@ void fft_3d(FFT_DATA *in, FFT_DATA *out, int flag, struct fft_plan_3d *plan)
total = plan->total3;
length = plan->length3;
#if defined(FFT_SGI)
for (offset = 0; offset < total; offset += length)
FFT_1D(flag,length,&data[offset],1,plan->coeff3);
#elif defined(FFT_SCSL)
for (offset = 0; offset < total; offset += length)
FFT_1D(flag,length,scalef,&data[offset],&data[offset],plan->coeff3,
plan->work3,&isys);
#elif defined(FFT_ACML)
num=total/length;
FFT_1D(&flag,&num,&length,data,plan->coeff3,&info);
#elif defined(FFT_INTEL)
for (offset = 0; offset < total; offset += length)
FFT_1D(&data[offset],&length,&flag,plan->coeff3);
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
if (flag == -1)
DftiComputeForward(plan->handle_slow,data);
else
DftiComputeBackward(plan->handle_slow,data);
#elif defined(FFT_DEC)
if (flag == -1)
for (offset = 0; offset < total; offset += length)
FFT_1D(&c,&c,&f,&data[offset],&data[offset],&length,&one);
else
for (offset = 0; offset < total; offset += length)
FFT_1D(&c,&c,&b,&data[offset],&data[offset],&length,&one);
#elif defined(FFT_T3E)
for (offset = 0; offset < total; offset += length)
FFT_1D(&flag,&length,&scalef,&data[offset],&data[offset],plan->coeff3,
plan->work3,&isys);
#elif defined(FFT_FFTW2)
if (flag == -1)
fftw(plan->plan_slow_forward,total/length,data,1,length,NULL,0,0);
@ -295,7 +210,6 @@ void fft_3d(FFT_DATA *in, FFT_DATA *out, int flag, struct fft_plan_3d *plan)
plan->post_plan);
// scaling if required
#if !defined(FFT_T3E) && !defined(FFT_ACML)
if (flag == 1 && plan->scaled) {
norm = plan->norm;
num = plan->normnum;
@ -314,25 +228,6 @@ void fft_3d(FFT_DATA *in, FFT_DATA *out, int flag, struct fft_plan_3d *plan)
#endif
}
}
#endif
#ifdef FFT_T3E
if (flag == 1 && plan->scaled) {
norm = plan->norm;
num = plan->normnum;
for (i = 0; i < num; i++) out[i] *= (norm,norm);
}
#endif
#ifdef FFT_ACML
norm = plan->norm;
num = plan->normnum;
for (i = 0; i < num; i++) {
out[i].re *= norm;
out[i].im *= norm;
}
#endif
}
/* ----------------------------------------------------------------------
@ -373,23 +268,6 @@ struct fft_plan_3d *fft_3d_create_plan(
int out_size,first_size,second_size,third_size,copy_size,scratch_size;
int np1,np2,ip1,ip2;
// system specific variables
#ifdef FFT_SCSL
FFT_DATA dummy_d[5];
FFT_PREC dummy_p[5];
int isign,isys;
FFT_PREC scalef;
#endif
#ifdef FFT_INTEL
FFT_DATA dummy;
#endif
#ifdef FFT_T3E
FFT_DATA dummy[5];
int isign,isys;
double scalef;
#endif
// query MPI info
MPI_Comm_rank(comm,&me);
@ -425,8 +303,7 @@ struct fft_plan_3d *fft_3d_create_plan(
first_klo = in_klo;
first_khi = in_khi;
plan->pre_plan = NULL;
}
else {
} else {
first_ilo = 0;
first_ihi = nfast - 1;
first_jlo = ip1*nmid/np1;
@ -487,8 +364,7 @@ struct fft_plan_3d *fft_3d_create_plan(
third_jhi = out_jhi;
third_klo = out_klo;
third_khi = out_khi;
}
else {
} else {
third_ilo = ip1*nfast/np1;
third_ihi = (ip1+1)*nfast/np1 - 1;
third_jlo = ip2*nmid/np2;
@ -605,139 +481,7 @@ struct fft_plan_3d *fft_3d_create_plan(
// system specific pre-computation of 1d FFT coeffs
// and scaling normalization
#if defined(FFT_SGI)
plan->coeff1 = (FFT_DATA *) malloc((nfast+15)*sizeof(FFT_DATA));
plan->coeff2 = (FFT_DATA *) malloc((nmid+15)*sizeof(FFT_DATA));
plan->coeff3 = (FFT_DATA *) malloc((nslow+15)*sizeof(FFT_DATA));
if (plan->coeff1 == NULL || plan->coeff2 == NULL ||
plan->coeff3 == NULL) return NULL;
FFT_1D_INIT(nfast,plan->coeff1);
FFT_1D_INIT(nmid,plan->coeff2);
FFT_1D_INIT(nslow,plan->coeff3);
if (scaled == 0)
plan->scaled = 0;
else {
plan->scaled = 1;
plan->norm = 1.0/(nfast*nmid*nslow);
plan->normnum = (out_ihi-out_ilo+1) * (out_jhi-out_jlo+1) *
(out_khi-out_klo+1);
}
#elif defined(FFT_SCSL)
plan->coeff1 = (FFT_PREC *) malloc((2*nfast+30)*sizeof(FFT_PREC));
plan->coeff2 = (FFT_PREC *) malloc((2*nmid+30)*sizeof(FFT_PREC));
plan->coeff3 = (FFT_PREC *) malloc((2*nslow+30)*sizeof(FFT_PREC));
if (plan->coeff1 == NULL || plan->coeff2 == NULL ||
plan->coeff3 == NULL) return NULL;
plan->work1 = (FFT_PREC *) malloc((2*nfast)*sizeof(FFT_PREC));
plan->work2 = (FFT_PREC *) malloc((2*nmid)*sizeof(FFT_PREC));
plan->work3 = (FFT_PREC *) malloc((2*nslow)*sizeof(FFT_PREC));
if (plan->work1 == NULL || plan->work2 == NULL ||
plan->work3 == NULL) return NULL;
isign = 0;
scalef = 1.0;
isys = 0;
FFT_1D_INIT(isign,nfast,scalef,dummy_d,dummy_d,plan->coeff1,dummy_p,&isys);
FFT_1D_INIT(isign,nmid,scalef,dummy_d,dummy_d,plan->coeff2,dummy_p,&isys);
FFT_1D_INIT(isign,nslow,scalef,dummy_d,dummy_d,plan->coeff3,dummy_p,&isys);
if (scaled == 0)
plan->scaled = 0;
else {
plan->scaled = 1;
plan->norm = 1.0/(nfast*nmid*nslow);
plan->normnum = (out_ihi-out_ilo+1) * (out_jhi-out_jlo+1) *
(out_khi-out_klo+1);
}
#elif defined(FFT_ACML)
plan->coeff1 = (FFT_DATA *) malloc((3*nfast+100)*sizeof(FFT_DATA));
plan->coeff2 = (FFT_DATA *) malloc((3*nmid+100)*sizeof(FFT_DATA));
plan->coeff3 = (FFT_DATA *) malloc((3*nslow+100)*sizeof(FFT_DATA));
if (plan->coeff1 == NULL || plan->coeff2 == NULL ||
plan->coeff3 == NULL) return NULL;
int isign = 100;
int isys = 1;
int info = 0;
FFT_DATA *dummy = NULL;
FFT_1D(&isign,&isys,&nfast,dummy,plan->coeff1,&info);
FFT_1D(&isign,&isys,&nmid,dummy,plan->coeff2,&info);
FFT_1D(&isign,&isys,&nslow,dummy,plan->coeff3,&info);
if (scaled == 0) {
plan->scaled = 0;
plan->norm = sqrt(nfast*nmid*nslow);
plan->normnum = (out_ihi-out_ilo+1) * (out_jhi-out_jlo+1) *
(out_khi-out_klo+1);
} else {
plan->scaled = 1;
plan->norm = sqrt(nfast*nmid*nslow);
plan->normnum = (out_ihi-out_ilo+1) * (out_jhi-out_jlo+1) *
(out_khi-out_klo+1);
}
#elif defined(FFT_INTEL)
flag = 0;
int i,num,fftflag;
int list[50];
num = 0;
factor(nfast,&num,list);
for (i = 0; i < num; i++)
if (list[i] != 2 && list[i] != 3 && list[i] != 5) flag = 1;
num = 0;
factor(nmid,&num,list);
for (i = 0; i < num; i++)
if (list[i] != 2 && list[i] != 3 && list[i] != 5) flag = 1;
num = 0;
factor(nslow,&num,list);
for (i = 0; i < num; i++)
if (list[i] != 2 && list[i] != 3 && list[i] != 5) flag = 1;
MPI_Allreduce(&flag,&fftflag,1,MPI_INT,MPI_MAX,comm);
if (fftflag) {
if (me == 0) printf("ERROR: FFTs are not power of 2,3,5\n");
return NULL;
}
plan->coeff1 = (FFT_DATA *) malloc((3*nfast/2+1)*sizeof(FFT_DATA));
plan->coeff2 = (FFT_DATA *) malloc((3*nmid/2+1)*sizeof(FFT_DATA));
plan->coeff3 = (FFT_DATA *) malloc((3*nslow/2+1)*sizeof(FFT_DATA));
if (plan->coeff1 == NULL || plan->coeff2 == NULL ||
plan->coeff3 == NULL) return NULL;
flag = 0;
FFT_1D_INIT(&dummy,&nfast,&flag,plan->coeff1);
FFT_1D_INIT(&dummy,&nmid,&flag,plan->coeff2);
FFT_1D_INIT(&dummy,&nslow,&flag,plan->coeff3);
if (scaled == 0) {
plan->scaled = 1;
plan->norm = nfast*nmid*nslow;
plan->normnum = (out_ihi-out_ilo+1) * (out_jhi-out_jlo+1) *
(out_khi-out_klo+1);
}
else
plan->scaled = 0;
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
DftiCreateDescriptor( &(plan->handle_fast), FFT_MKL_PREC, DFTI_COMPLEX, 1, (MKL_LONG)nfast);
DftiSetValue(plan->handle_fast, DFTI_NUMBER_OF_TRANSFORMS, (MKL_LONG)plan->total1/nfast);
DftiSetValue(plan->handle_fast, DFTI_PLACEMENT,DFTI_INPLACE);
@ -768,50 +512,6 @@ struct fft_plan_3d *fft_3d_create_plan(
(out_khi-out_klo+1);
}
#elif defined(FFT_DEC)
if (scaled == 0) {
plan->scaled = 1;
plan->norm = nfast*nmid*nslow;
plan->normnum = (out_ihi-out_ilo+1) * (out_jhi-out_jlo+1) *
(out_khi-out_klo+1);
}
else
plan->scaled = 0;
#elif defined(FFT_T3E)
plan->coeff1 = (double *) malloc((12*nfast)*sizeof(double));
plan->coeff2 = (double *) malloc((12*nmid)*sizeof(double));
plan->coeff3 = (double *) malloc((12*nslow)*sizeof(double));
if (plan->coeff1 == NULL || plan->coeff2 == NULL ||
plan->coeff3 == NULL) return NULL;
plan->work1 = (double *) malloc((8*nfast)*sizeof(double));
plan->work2 = (double *) malloc((8*nmid)*sizeof(double));
plan->work3 = (double *) malloc((8*nslow)*sizeof(double));
if (plan->work1 == NULL || plan->work2 == NULL ||
plan->work3 == NULL) return NULL;
isign = 0;
scalef = 1.0;
isys = 0;
FFT_1D_INIT(&isign,&nfast,&scalef,dummy,dummy,plan->coeff1,dummy,&isys);
FFT_1D_INIT(&isign,&nmid,&scalef,dummy,dummy,plan->coeff2,dummy,&isys);
FFT_1D_INIT(&isign,&nslow,&scalef,dummy,dummy,plan->coeff3,dummy,&isys);
if (scaled == 0)
plan->scaled = 0;
else {
plan->scaled = 1;
plan->norm = 1.0/(nfast*nmid*nslow);
plan->normnum = (out_ihi-out_ilo+1) * (out_jhi-out_jlo+1) *
(out_khi-out_klo+1);
}
#elif defined(FFT_FFTW2)
plan->plan_fast_forward =
@ -948,36 +648,10 @@ void fft_3d_destroy_plan(struct fft_plan_3d *plan)
if (plan->copy) free(plan->copy);
if (plan->scratch) free(plan->scratch);
#if defined(FFT_SGI)
free(plan->coeff1);
free(plan->coeff2);
free(plan->coeff3);
#elif defined(FFT_SCSL)
free(plan->coeff1);
free(plan->coeff2);
free(plan->coeff3);
free(plan->work1);
free(plan->work2);
free(plan->work3);
#elif defined(FFT_ACML)
free(plan->coeff1);
free(plan->coeff2);
free(plan->coeff3);
#elif defined(FFT_INTEL)
free(plan->coeff1);
free(plan->coeff2);
free(plan->coeff3);
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
DftiFreeDescriptor(&(plan->handle_fast));
DftiFreeDescriptor(&(plan->handle_mid));
DftiFreeDescriptor(&(plan->handle_slow));
#elif defined(FFT_T3E)
free(plan->coeff1);
free(plan->coeff2);
free(plan->coeff3);
free(plan->work1);
free(plan->work2);
free(plan->work3);
#elif defined(FFT_FFTW2)
if (plan->plan_slow_forward != plan->plan_fast_forward &&
plan->plan_slow_forward != plan->plan_mid_forward) {
@ -1022,38 +696,31 @@ void factor(int n, int *num, int *list)
{
if (n == 1) {
return;
}
else if (n % 2 == 0) {
} else if (n % 2 == 0) {
*list = 2;
(*num)++;
factor(n/2,num,list+1);
}
else if (n % 3 == 0) {
} else if (n % 3 == 0) {
*list = 3;
(*num)++;
factor(n/3,num,list+1);
}
else if (n % 5 == 0) {
} else if (n % 5 == 0) {
*list = 5;
(*num)++;
factor(n/5,num,list+1);
}
else if (n % 7 == 0) {
} else if (n % 7 == 0) {
*list = 7;
(*num)++;
factor(n/7,num,list+1);
}
else if (n % 11 == 0) {
} else if (n % 11 == 0) {
*list = 11;
(*num)++;
factor(n/11,num,list+1);
}
else if (n % 13 == 0) {
} else if (n % 13 == 0) {
*list = 13;
(*num)++;
factor(n/13,num,list+1);
}
else {
} else {
*list = n;
(*num)++;
return;
@ -1099,23 +766,6 @@ void fft_1d_only(FFT_DATA *data, int nsize, int flag, struct fft_plan_3d *plan)
FFT_SCALAR *data_ptr;
#endif
// system specific constants
#ifdef FFT_SCSL
int isys = 0;
FFT_PREC scalef = 1.0;
#endif
#ifdef FFT_DEC
char c = 'C';
char f = 'F';
char b = 'B';
int one = 1;
#endif
#ifdef FFT_T3E
int isys = 0;
double scalef = 1.0;
#endif
// total = size of data needed in each dim
// length = length of 1d FFT in each dim
// total/length = # of 1d FFTs in each dim
@ -1141,39 +791,7 @@ void fft_1d_only(FFT_DATA *data, int nsize, int flag, struct fft_plan_3d *plan)
// perform 1d FFTs in each of 3 dimensions
// data is just an array of 0.0
#ifdef FFT_SGI
for (int offset = 0; offset < total1; offset += length1)
FFT_1D(flag,length1,&data[offset],1,plan->coeff1);
for (int offset = 0; offset < total2; offset += length2)
FFT_1D(flag,length2,&data[offset],1,plan->coeff2);
for (int offset = 0; offset < total3; offset += length3)
FFT_1D(flag,length3,&data[offset],1,plan->coeff3);
#elif defined(FFT_SCSL)
for (int offset = 0; offset < total1; offset += length1)
FFT_1D(flag,length1,scalef,&data[offset],&data[offset],plan->coeff1,
plan->work1,&isys);
for (int offset = 0; offset < total2; offset += length2)
FFT_1D(flag,length2,scalef,&data[offset],&data[offset],plan->coeff2,
plan->work2,&isys);
for (int offset = 0; offset < total3; offset += length3)
FFT_1D(flag,length3,scalef,&data[offset],&data[offset],plan->coeff3,
plan->work3,&isys);
#elif defined(FFT_ACML)
int info=0;
num=total1/length1;
FFT_1D(&flag,&num,&length1,data,plan->coeff1,&info);
num=total2/length2;
FFT_1D(&flag,&num,&length2,data,plan->coeff2,&info);
num=total3/length3;
FFT_1D(&flag,&num,&length3,data,plan->coeff3,&info);
#elif defined(FFT_INTEL)
for (int offset = 0; offset < total1; offset += length1)
FFT_1D(&data[offset],&length1,&flag,plan->coeff1);
for (int offset = 0; offset < total2; offset += length2)
FFT_1D(&data[offset],&length2,&flag,plan->coeff2);
for (int offset = 0; offset < total3; offset += length3)
FFT_1D(&data[offset],&length3,&flag,plan->coeff3);
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
if (flag == -1) {
DftiComputeForward(plan->handle_fast,data);
DftiComputeForward(plan->handle_mid,data);
@ -1183,32 +801,6 @@ void fft_1d_only(FFT_DATA *data, int nsize, int flag, struct fft_plan_3d *plan)
DftiComputeBackward(plan->handle_mid,data);
DftiComputeBackward(plan->handle_slow,data);
}
#elif defined(FFT_DEC)
if (flag == -1) {
for (int offset = 0; offset < total1; offset += length1)
FFT_1D(&c,&c,&f,&data[offset],&data[offset],&length1,&one);
for (int offset = 0; offset < total2; offset += length2)
FFT_1D(&c,&c,&f,&data[offset],&data[offset],&length2,&one);
for (int offset = 0; offset < total3; offset += length3)
FFT_1D(&c,&c,&f,&data[offset],&data[offset],&length3,&one);
} else {
for (int offset = 0; offset < total1; offset += length1)
FFT_1D(&c,&c,&b,&data[offset],&data[offset],&length1,&one);
for (int offset = 0; offset < total2; offset += length2)
FFT_1D(&c,&c,&b,&data[offset],&data[offset],&length2,&one);
for (int offset = 0; offset < total3; offset += length3)
FFT_1D(&c,&c,&b,&data[offset],&data[offset],&length3,&one);
}
#elif defined(FFT_T3E)
for (int offset = 0; offset < total1; offset += length1)
FFT_1D(&flag,&length1,&scalef,&data[offset],&data[offset],plan->coeff1,
plan->work1,&isys);
for (int offset = 0; offset < total2; offset += length2)
FFT_1D(&flag,&length2,&scalef,&data[offset],&data[offset],plan->coeff2,
plan->work2,&isys);
for (int offset = 0; offset < total3; offset += length3)
FFT_1D(&flag,&length3,&scalef,&data[offset],&data[offset],plan->coeff3,
plan->work3,&isys);
#elif defined(FFT_FFTW2)
if (flag == -1) {
fftw(plan->plan_fast_forward,total1/length1,data,1,0,NULL,0,0);
@ -1257,7 +849,6 @@ void fft_1d_only(FFT_DATA *data, int nsize, int flag, struct fft_plan_3d *plan)
// scaling if required
// limit num to size of data
#ifndef FFT_T3E
if (flag == 1 && plan->scaled) {
norm = plan->norm;
num = MIN(plan->normnum,nsize);
@ -1276,13 +867,4 @@ void fft_1d_only(FFT_DATA *data, int nsize, int flag, struct fft_plan_3d *plan)
#endif
}
}
#endif
#ifdef FFT_T3E
if (flag == 1 && plan->scaled) {
norm = plan->norm;
num = MIN(plan->normnum,nsize);
for (i = 0; i < num; i++) data[i] *= (norm,norm);
}
#endif
}

158
src/KSPACE/fft3d.h
View File

@ -27,7 +27,7 @@ typedef double FFT_SCALAR;
// set default fftw library. switch to FFT_FFTW3 when convenient.
#ifdef FFT_FFTW
#define FFT_FFTW2
#define FFT_FFTW3
#endif
// -------------------------------------------------------------------------
@ -36,73 +36,11 @@ typedef double FFT_SCALAR;
#if FFT_PRECISION == 1
#if defined(FFT_SGI)
#include "fft.h"
typedef complex FFT_DATA;
#define FFT_1D cfft1d
#define FFT_1D_INIT cfft1di
extern "C" {
int cfft1d(int, int, FFT_DATA *, int, FFT_DATA *);
FFT_DATA *cfft1di(int, FFT_DATA *);
}
#elif defined(FFT_SCSL)
#include <scsl_fft.h>
typedef scsl_complex FFT_DATA;
typedef float FFT_PREC;
#define FFT_1D ccfft
#define FFT_1D_INIT ccfft
extern "C" {
int ccfft(int, int, FFT_PREC, FFT_DATA *, FFT_DATA *,
FFT_PREC *, FFT_PREC *, int *);
}
#elif defined(FFT_ACML)
typedef struct {
float re;
float im;
} FFT_DATA;
#define FFT_1D cfft1m_
extern "C" {
void cfft1m_(int *, int *, int *, FFT_DATA *, FFT_DATA *, int *);
}
#elif defined(FFT_INTEL)
typedef struct {
float re;
float im;
} FFT_DATA;
#define FFT_1D cfft1d_
#define FFT_1D_INIT cfft1d_
extern "C" {
void cfft1d_(FFT_DATA *, int *, int *, FFT_DATA *);
}
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
#include "mkl_dfti.h"
typedef float _Complex FFT_DATA;
#define FFT_MKL_PREC DFTI_SINGLE
#elif defined(FFT_DEC)
typedef struct {
float re;
float im;
} FFT_DATA;
#define FFT_1D cfft_
extern "C" {
void cfft_(char *, char *, char *, FFT_DATA *, FFT_DATA *, int *, int *);
}
#elif defined(FFT_T3E)
#include <complex.h>
typedef complex single FFT_DATA;
#define FFT_1D GGFFT
#define FFT_1D_INIT GGFFT
extern "C" {
void GGFFT(int *, int *, double *, FFT_DATA *, FFT_DATA *,
double *, double *, int *);
}
#elif defined(FFT_FFTW2)
#if defined(FFTW_SIZE)
#include "sfftw.h"
@ -138,73 +76,11 @@ typedef struct kiss_fft_state* kiss_fft_cfg;
#elif FFT_PRECISION == 2
#if defined(FFT_SGI)
#include "fft.h"
typedef zomplex FFT_DATA;
#define FFT_1D zfft1d
#define FFT_1D_INIT zfft1di
extern "C" {
int zfft1d(int, int, FFT_DATA *, int, FFT_DATA *);
FFT_DATA *zfft1di(int, FFT_DATA *);
}
#elif defined(FFT_SCSL)
#include <scsl_fft.h>
typedef scsl_zomplex FFT_DATA;
typedef double FFT_PREC;
#define FFT_1D zzfft
#define FFT_1D_INIT zzfft
extern "C" {
int zzfft(int, int, FFT_PREC, FFT_DATA *, FFT_DATA *,
FFT_PREC *, FFT_PREC *, int *);
}
#elif defined(FFT_ACML)
typedef struct {
double re;
double im;
} FFT_DATA;
#define FFT_1D zfft1m_
extern "C" {
void zfft1m_(int *, int *, int *, FFT_DATA *, FFT_DATA *, int *);
}
#elif defined(FFT_INTEL)
typedef struct {
double re;
double im;
} FFT_DATA;
#define FFT_1D zfft1d_
#define FFT_1D_INIT zfft1d_
extern "C" {
void zfft1d_(FFT_DATA *, int *, int *, FFT_DATA *);
}
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
#include "mkl_dfti.h"
typedef double _Complex FFT_DATA;
#define FFT_MKL_PREC DFTI_DOUBLE
#elif defined(FFT_DEC)
typedef struct {
double re;
double im;
} FFT_DATA;
#define FFT_1D zfft_
extern "C" {
void zfft_(char *, char *, char *, FFT_DATA *, FFT_DATA *, int *, int *);
}
#elif defined(FFT_T3E)
#include <complex.h>
typedef complex double FFT_DATA;
#define FFT_1D CCFFT
#define FFT_1D_INIT CCFFT
extern "C" {
void CCFFT(int *, int *, double *, FFT_DATA *, FFT_DATA *,
double *, double *, int *);
}
#elif defined(FFT_FFTW2)
#if defined(FFTW_SIZE)
#include "dfftw.h"
@ -258,36 +134,10 @@ struct fft_plan_3d {
double norm; // normalization factor for rescaling
// system specific 1d FFT info
#if defined(FFT_SGI)
FFT_DATA *coeff1;
FFT_DATA *coeff2;
FFT_DATA *coeff3;
#elif defined(FFT_SCSL)
FFT_PREC *coeff1;
FFT_PREC *coeff2;
FFT_PREC *coeff3;
FFT_PREC *work1;
FFT_PREC *work2;
FFT_PREC *work3;
#elif defined(FFT_ACML)
FFT_DATA *coeff1;
FFT_DATA *coeff2;
FFT_DATA *coeff3;
#elif defined(FFT_INTEL)
FFT_DATA *coeff1;
FFT_DATA *coeff2;
FFT_DATA *coeff3;
#elif defined(FFT_MKL)
#if defined(FFT_MKL)
DFTI_DESCRIPTOR *handle_fast;
DFTI_DESCRIPTOR *handle_mid;
DFTI_DESCRIPTOR *handle_slow;
#elif defined(FFT_T3E)
double *coeff1;
double *coeff2;
double *coeff3;
double *work1;
double *work2;
double *work3;
#elif defined(FFT_FFTW2)
fftw_plan plan_fast_forward;
fftw_plan plan_fast_backward;

2
src/MAKE/MACHINES/Makefile.bgl
View File

@ -56,7 +56,7 @@ MPI_LIB = -lmpich.rts
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW
FFT_INC = -DFFT_FFTW2
FFT_PATH =
FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.bgq
View File

@ -8,7 +8,7 @@ SHELL = /bin/bash
# do not edit this section
# select which compiler by editing Makefile.bgq.details
include ../MAKE/Makefile.bgq.details
include ../MAKE/MACHINES/Makefile.bgq.details
include Makefile.package.settings
include Makefile.package

2
src/MAKE/MACHINES/Makefile.bgq.details
View File

@ -101,7 +101,7 @@ MPI_LIB += #/home/jhammond/OSPRI/branches/marpn/wrap/libmpiarbrpn.a
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -I/soft/libraries/alcf/current/xl/FFTW2/include -DFFT_FFTW -DFFTW_SIZE
FFT_INC = -I/soft/libraries/alcf/current/xl/FFTW2/include -DFFT_FFTW2 -DFFTW_SIZE
FFT_PATH = #/soft/libraries/alcf/current/xl/FFTW2
FFT_LIB = -L/soft/libraries/alcf/current/xl/FFTW2/lib -ldfftw

2
src/MAKE/MACHINES/Makefile.chama
View File

@ -50,7 +50,7 @@ MPI_LIB =
# PATH = path for FFT library
# LIB = name of FFT library
#FFT_INC = -DFFT_FFTW -I${FFTW_INCLUDE}
#FFT_INC = -DFFT_FFTW2 -I${FFTW_INCLUDE}
#FFT_PATH = -L${FFTW_LIB}
#FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.cygwin
View File

@ -50,7 +50,7 @@ MPI_LIB = -lmpich
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW -I/cygdrive/c/cygwin/usr/local/include
FFT_INC = -DFFT_FFTW2 -I/cygdrive/c/cygwin/usr/local/include
FFT_PATH = -L/cygdrive/c/cygwin/usr/local/lib
FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.glory
View File

@ -67,7 +67,7 @@ MPI_LIB =
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW -I${FFTW_INCLUDE}
FFT_INC = -DFFT_FFTW2 -I${FFTW_INCLUDE}
FFT_PATH = -L${FFTW_LIB}
FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.mac
View File

@ -50,7 +50,7 @@ MPI_LIB = -lmpi_stubs
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW
FFT_INC = -DFFT_FFTW2
FFT_PATH =
FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.mac_mpi
View File

@ -53,7 +53,7 @@ MPI_LIB =
FFTW = /usr/local
FFT_INC = -DFFT_FFTW -I${FFTW}/include
FFT_INC = -DFFT_FFTW2 -I${FFTW}/include
FFT_PATH = -L${FFTW}/lib
FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.power
View File

@ -51,7 +51,7 @@ MPI_LIB =
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW -I/scr/oppe/LAMMPS/fftw-2.1.5/include
FFT_INC = -DFFT_FFTW2 -I/scr/oppe/LAMMPS/fftw-2.1.5/include
FFT_PATH = -L/scr/oppe/LAMMPS/fftw-2.1.5/lib
FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.redsky
View File

@ -70,7 +70,7 @@ MPI_LIB =
# PATH = path for FFT library
# LIB = name of FFT library
#FFT_INC = -DFFT_FFTW -I${FFTW_INCLUDE}
#FFT_INC = -DFFT_FFTW2 -I${FFTW_INCLUDE}
#FFT_PATH = -L${FFTW_LIB}
#FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.storm
View File

@ -51,7 +51,7 @@ MPI_LIB =
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW -I/projects/fftw/fftw-2.1.5/include
FFT_INC = -DFFT_FFTW2 -I/projects/fftw/fftw-2.1.5/include
FFT_PATH =
FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.tacc
View File

@ -53,7 +53,7 @@ MPI_LIB = -lmpich -lpthread
FFTW_INC = ${TACC_FFTW2_INC}
FFTW_LIB = ${TACC_FFTW2_LIB}
FFT_INC = -DFFT_FFTW -I${FFTW_INC}
FFT_INC = -DFFT_FFTW2 -I${FFTW_INC}
FFT_PATH = -L${FFTW_LIB}
FFT_LIB = ${FFTW_LIB}/libfftw.a

2
src/MAKE/MACHINES/Makefile.xe6
View File

@ -51,7 +51,7 @@ MPI_LIB =
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW
FFT_INC = -DFFT_FFTW2
FFT_PATH =
FFT_LIB = -lfftw

2
src/MAKE/MACHINES/Makefile.xt3
View File

@ -52,7 +52,7 @@ MPI_LIB = -lmpich -lpthread
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW -I$(FFTW_INC)
FFT_INC = -DFFT_FFTW2 -I$(FFTW_INC)
FFT_PATH = -L$(FFTW_LIB)
FFT_LIB = -ldfftw

2
src/MAKE/MACHINES/Makefile.xt5
View File

@ -51,7 +51,7 @@ MPI_LIB =
# PATH = path for FFT library
# LIB = name of FFT library
FFT_INC = -DFFT_FFTW -I/home/sjplimp/fftw/fftw
FFT_INC = -DFFT_FFTW2 -I/home/sjplimp/fftw/fftw
FFT_PATH = -L/home/sjplimp/fftw/fftw/.libs
FFT_LIB = -lfftw