Commit e70a251fdd863e281bbf3ae30319e4ed8ae95a97

Authored by dmayerich
1 parent 20e5e34a

fixed double precision errors

Showing 3 changed files with 35 additions and 6 deletions   Show diff stats
fieldslice.cpp
... ... @@ -23,19 +23,27 @@ fieldslice::fieldslice(unsigned int x_size, unsigned int y_size)
23 23 void fieldslice::toAngularSpectrum()
24 24 {
25 25 cufftHandle plan;
  26 + cufftResult result;
26 27  
27 28 //create a CUFFT plan handle
28   - if(cufftPlan2d(&plan, R[0], R[1], CUFFT_C2C) != CUFFT_SUCCESS)
  29 +#ifdef PRECISION_SINGLE
  30 + result = cufftPlan2d(&plan, R[0], R[1], CUFFT_C2C);
  31 +#elif defined PRECISION_DOUBLE
  32 + result = cufftPlan2d(&plan, R[0], R[1], CUFFT_Z2Z);
  33 +#endif
  34 +
  35 + if(result != CUFFT_SUCCESS)
29 36 {
30 37 cout<<"Error creating CUFFT plan for computing the angular spectrum."<<endl;
31 38 exit(1);
32 39 }
33 40  
34 41 #ifdef PRECISION_SINGLE
35   - if(cufftExecC2C(plan, (cufftComplex*)x_hat, (cufftComplex*)x_hat, CUFFT_FORWARD) != CUFFT_SUCCESS)
  42 + result = cufftExecC2C(plan, (cufftComplex*)x_hat, (cufftComplex*)x_hat, CUFFT_FORWARD);
36 43 #elif defined PRECISION_DOUBLE
37   - if(cufftExecZ2Z(plan, (cufftDoubleComplex*)x_hat, (cufftDoubleComplex*)x_hat, CUFFT_FORWARD) != CUFFT_SUCCESS)
  44 + result = cufftExecZ2Z(plan, (cufftDoubleComplex*)x_hat, (cufftDoubleComplex*)x_hat, CUFFT_FORWARD);
38 45 #endif
  46 + if(result != CUFFT_SUCCESS)
39 47 {
40 48 cout<<"Error executing the CUFFT forward FFT to compute the angular spectrum."<<endl;
41 49 exit(1);
... ... @@ -49,19 +57,27 @@ void fieldslice::toAngularSpectrum()
49 57 void fieldslice::fromAngularSpectrum()
50 58 {
51 59 cufftHandle plan;
  60 + cufftResult result;
52 61  
53 62 //create a CUFFT plan handle
54   - if(cufftPlan2d(&plan, R[0], R[1], CUFFT_C2C) != CUFFT_SUCCESS)
  63 +#ifdef PRECISION_SINGLE
  64 + result = cufftPlan2d(&plan, R[0], R[1], CUFFT_C2C);
  65 +#elif defined PRECISION_DOUBLE
  66 + result = cufftPlan2d(&plan, R[0], R[1], CUFFT_Z2Z);
  67 +#endif
  68 +
  69 + if(result != CUFFT_SUCCESS)
55 70 {
56 71 cout<<"Error creating CUFFT plan for computing the angular spectrum."<<endl;
57 72 exit(1);
58 73 }
59 74  
60 75 #ifdef PRECISION_SINGLE
61   - if(cufftExecC2C(plan, (cufftComplex*)x_hat, (cufftComplex*)x_hat, CUFFT_INVERSE) != CUFFT_SUCCESS)
  76 + result = cufftExecC2C(plan, (cufftComplex*)x_hat, (cufftComplex*)x_hat, CUFFT_INVERSE);
62 77 #elif defined PRECISION_DOUBLE
63   - if(cufftExecZ2Z(plan, (cufftDoubleComplex*)x_hat, (cufftDoubleComplex*)x_hat, CUFFT_INVERSE) != CUFFT_SUCCESS)
  78 + result = cufftExecZ2Z(plan, (cufftDoubleComplex*)x_hat, (cufftDoubleComplex*)x_hat, CUFFT_INVERSE);
64 79 #endif
  80 + if(result != CUFFT_SUCCESS)
65 81 {
66 82 cout<<"Error executing the CUFFT forward FFT to compute the angular spectrum."<<endl;
67 83 exit(1);
... ...
scalarslice.cu
... ... @@ -47,7 +47,20 @@ void scalarslice::toImage(std::string filename, bool positive, rts::colormap::co
47 47 //find the index of the value with maximum magnitude
48 48 int N = R[0] * R[1];
49 49 int result;
  50 +#ifdef PRECISION_SINGLE
50 51 stat = cublasIsamax(handle, N, S, 1, &result);
  52 +#elif defined PRECISION_DOUBLE
  53 + stat = cublasIdamax(handle, N, S, 1, &result);
  54 +#endif
  55 +
  56 + //adjust for 1-based indexing
  57 + result -= 1;
  58 +
  59 + if(stat != CUBLAS_STATUS_SUCCESS)
  60 + {
  61 + std::cout<<"CUBLAS Error: failure finding maximum value."<<std::endl;
  62 + exit(1);
  63 + }
51 64  
52 65 //std::cout<<"Maximum index: "<<result<<std::endl;
53 66  
... ...
win32/bimsim.exe
No preview for this file type