![]() |
kspaceFirstOrder3D-OMP 1.0
The C++ implementation of the k-wave toolbox for the time-domain simulation of acoustic wave fields in 3D
|
00001 00032 #include <MatrixClasses/FFTWComplexMatrix.h> 00033 #include <MatrixClasses/RealMatrix.h> 00034 00035 #include <Utils/ErrorMessages.h> 00036 00037 using namespace std; 00038 00039 00040 00041 //----------------------------------------------------------------------------// 00042 // Constants // 00043 //----------------------------------------------------------------------------// 00044 00045 //----------------------------------------------------------------------------// 00046 // Public methods // 00047 //----------------------------------------------------------------------------// 00048 00049 00054 TFFTWComplexMatrix::TFFTWComplexMatrix(struct TDimensionSizes DimensionSizes) : 00055 TComplexMatrix(), 00056 fftw_plan_R2C_Created(false), fftw_plan_C2R_Created(false) 00057 { 00058 00059 InitDimensions(DimensionSizes); 00060 00061 AllocateMemory(); 00062 }// end of TFFTWComplexMatrix 00063 //------------------------------------------------------------------------------ 00064 00065 00066 00070 TFFTWComplexMatrix::~TFFTWComplexMatrix(){ 00071 00072 //-- Destroy fftw plans --// 00073 00074 FreeMemory(); 00075 00076 if (fftw_plan_R2C_Created) fftwf_destroy_plan(fftw_plan_R2C); 00077 if (fftw_plan_C2R_Created) fftwf_destroy_plan(fftw_plan_C2R); 00078 00079 00080 }// end of ~TFFTWComplexMatrix() 00081 //------------------------------------------------------------------------------ 00082 00083 00084 00090 void TFFTWComplexMatrix::CreateFFTPlan3D_R2C(TRealMatrix& InMatrix){ 00091 00092 00093 fftw_plan_R2C = fftwf_plan_dft_r2c_3d(InMatrix.GetDimensionSizes().Z,InMatrix.GetDimensionSizes().Y,InMatrix.GetDimensionSizes().X, 00094 InMatrix.GetRawData(), (fftwf_complex *) pMatrixData, 00095 TFFTWComplexMatrix_FFT_FLAG); 00096 00097 fftw_plan_R2C_Created = true; 00098 00099 00100 }// end of CreateFFTPlan3D_RealToComplex 00101 //------------------------------------------------------------------------------ 00102 00108 void TFFTWComplexMatrix::CreateFFTPlan3D_C2R(TRealMatrix& OutMatrix){ 00109 00110 00111 fftw_plan_C2R = fftwf_plan_dft_c2r_3d(OutMatrix.GetDimensionSizes().Z,OutMatrix.GetDimensionSizes().Y,OutMatrix.GetDimensionSizes().X, 00112 (fftwf_complex *)pMatrixData, OutMatrix.GetRawData(), 00113 TFFTWComplexMatrix_FFT_FLAG); 00114 00115 fftw_plan_C2R_Created = true; 00116 00117 }//end of CreateFFTPlan3D_ComplexToReal 00118 //------------------------------------------------------------------------------ 00119 00120 00121 00122 00127 void TFFTWComplexMatrix::Compute_FFT_3D_R2C(TRealMatrix & InMatrix){ 00128 00129 00130 // execute real to complex fft 00131 fftwf_execute_dft_r2c(fftw_plan_R2C, InMatrix.GetRawData(), (fftwf_complex *) pMatrixData); 00132 00133 00134 00135 }// end of Compute_FFT_3D_r2c 00136 //------------------------------------------------------------------------------ 00137 00142 void TFFTWComplexMatrix::Compute_iFFT_3D_C2R(TRealMatrix & OutMatrix){ 00143 00144 fftwf_execute_dft_c2r(fftw_plan_C2R,(fftwf_complex *) pMatrixData, OutMatrix.GetRawData()); 00145 00146 }// end of Compute_iFFT_3D_c2r 00147 //------------------------------------------------------------------------------ 00148 00149 00150 00151 //----------------------------------------------------------------------------// 00152 // Protected methods // 00153 //----------------------------------------------------------------------------// 00154 00159 void TFFTWComplexMatrix::AllocateMemory(){ 00160 00161 /* No memory allocated before this function*/ 00162 00163 00164 pMatrixData = (float *) fftwf_malloc(pTotalAllocatedElementCount * sizeof (float)); 00165 00166 00167 if (!pMatrixData) { 00168 fprintf(stderr,Matrix_ERR_FMT_NotEnoughMemory, "TFFTWComplexMatrix"); 00169 throw bad_alloc(); 00170 00171 } 00172 00173 //-- first touch --// 00174 #ifndef __NO_OMP__ 00175 #pragma omp parallel for schedule(static) 00176 #endif 00177 for (size_t i=0; i<pTotalAllocatedElementCount; i++){ 00178 pMatrixData[i] = 0.0f; 00179 } 00180 00181 00182 }// end of virtual bool FreeMemory 00183 //----------------------------------------------------------------------------- 00184 00185 00189 void TFFTWComplexMatrix::FreeMemory(){ 00190 00191 if (pMatrixData) fftwf_free( pMatrixData); 00192 pMatrixData = NULL; 00193 00194 }// end of FreeMemory 00195 //----------------------------------------------------------------------------- 00196 00197 00198 00199 00200 00201 00202 //----------------------------------------------------------------------------// 00203 // Private methods // 00204 //----------------------------------------------------------------------------//