![]() |
kspaceFirstOrder3D-OMP 1.0
The C++ implementation of the k-wave toolbox for the time-domain simulation of acoustic wave fields in 3D
|
00001 00031 #ifndef MATRIXCONTAINER_H 00032 #define MATRIXCONTAINER_H 00033 00034 #include <string.h> 00035 #include <map> 00036 00037 #include <MatrixClasses/BaseMatrix.h> 00038 #include <MatrixClasses/BaseFloatMatrix.h> 00039 #include <MatrixClasses/RealMatrix.h> 00040 #include <MatrixClasses/ComplexMatrix.h> 00041 #include <MatrixClasses/FFTWComplexMatrix.h> 00042 #include <MatrixClasses/UXYZ_SGXYZMatrix.h> 00043 #include <MatrixClasses/LongMatrix.h> 00044 00045 #include <Utils/MatrixNames.h> 00046 #include <Utils/DimensionSizes.h> 00047 00048 00053 enum TMatrixID {kappa, c2, p, 00054 00055 ux_sgx,uy_sgy, uz_sgz, 00056 duxdx, duydy, duzdz, 00057 dxudxn , dyudyn , dzudzn, 00058 dxudxn_sgx, dyudyn_sgy, dzudzn_sgz, 00059 00060 rhox, rhoy, rhoz, rho0, 00061 dt_rho0_sgx, dt_rho0_sgy, dt_rho0_sgz, 00062 00063 00064 p0_source_input, sensor_mask_ind, 00065 ddx_k_shift_pos, ddy_k_shift_pos, ddz_k_shift_pos, 00066 ddx_k_shift_neg, ddy_k_shift_neg, ddz_k_shift_neg, 00067 00068 pml_x_sgx, pml_y_sgy, pml_z_sgz, 00069 pml_x , pml_y , pml_z, 00070 00071 absorb_tau, absorb_eta, absorb_nabla1, absorb_nabla2, BonA, 00072 00073 ux_source_input, uy_source_input, uz_source_input, 00074 p_source_input, 00075 00076 u_source_index, p_source_index, transducer_source_input, 00077 delay_mask, 00078 00079 //---------------- output matrices -------------// 00080 p_sensor_raw, p_sensor_rms, p_sensor_max, p_sensor_i_1_raw, 00081 ux_sensor_raw, uy_sensor_raw, uz_sensor_raw, 00082 ux_sensor_i_1_agr_2, uy_sensor_i_1_agr_2, uz_sensor_i_1_agr_2,// aggregated values form two points 00083 00084 ux_sensor_rms, uy_sensor_rms, uz_sensor_rms, 00085 ux_sensor_max, uy_sensor_max, uz_sensor_max, 00086 00087 Ix_sensor_avg, Iy_sensor_avg, Iz_sensor_avg, 00088 Ix_sensor_max, Iy_sensor_max, Iz_sensor_max, 00089 00090 00091 //--------------Temporary matrices -------------// 00092 Temp_1_RS3D, Temp_2_RS3D, Temp_3_RS3D, 00093 FFT_X_temp, FFT_Y_temp, FFT_Z_temp 00094 }; 00095 00096 00097 00103 struct TMatrixRecord{ 00104 00109 enum TMatrixDataType {mdtReal, mdtComplex, mdtIndex, mdtFFTW, mdtUxyz}; 00110 00112 TBaseMatrix * MatrixPtr; 00114 TMatrixDataType MatrixDataType; 00116 TDimensionSizes DimensionSizes; 00118 bool LoadData; 00120 string HDF5MatrixName; 00121 00123 TMatrixRecord() : MatrixPtr(NULL), MatrixDataType(mdtReal), 00124 DimensionSizes(), LoadData(false), HDF5MatrixName("") 00125 {}; 00126 00128 TMatrixRecord(const TMatrixRecord& src); 00129 00131 TMatrixRecord& operator = (const TMatrixRecord& src); 00132 00134 void SetAllValues(TBaseMatrix * MatrixPtr, 00135 const TMatrixDataType MatrixDataType, 00136 const TDimensionSizes DimensionSizes, 00137 const bool LoadData, 00138 const string HDF5MatrixName); 00139 00140 virtual ~TMatrixRecord() {}; 00141 00142 00143 };// end of TMatrixRecord 00144 //------------------------------------------------------------------------------ 00145 00146 00147 00152 class TMatrixContainer { 00153 public: 00154 00156 TMatrixContainer() {} 00158 virtual ~TMatrixContainer(); 00159 00164 size_t size () {return MatrixContainer.size(); }; 00169 bool empty() {return MatrixContainer.empty();}; 00170 00171 00173 void CreateAllObjects(); 00175 void LoadMatricesDataFromDisk(THDF5_File & HDF5_File); 00177 void FreeAllMatrices(); 00178 00179 00181 void AddMatricesIntoContainer(); 00182 00188 TMatrixRecord& GetMatrixRecord (const TMatrixID MatrixID){ 00189 return MatrixContainer[MatrixID]; 00190 }; 00191 00197 TMatrixRecord& operator [] (const TMatrixID MatrixID){ 00198 return MatrixContainer[MatrixID]; 00199 }; 00200 00201 00207 TBaseMatrix & GetBaseMatrix (const TMatrixID MatrixID){ 00208 return static_cast<TBaseMatrix &> (*(MatrixContainer[MatrixID].MatrixPtr)); 00209 }; 00210 00216 TBaseFloatMatrix & GetBaseFloatMatrix (const TMatrixID MatrixID){ 00217 return static_cast<TBaseFloatMatrix &> (*(MatrixContainer[MatrixID].MatrixPtr)); 00218 }; 00219 00225 TRealMatrix & GetRealMatrix (const TMatrixID MatrixID){ 00226 return static_cast<TRealMatrix &> (*(MatrixContainer[MatrixID].MatrixPtr)); 00227 }; 00228 00234 Tuxyz_sgxyzMatrix & GetUxyz_sgxyzMatrix (const TMatrixID MatrixID){ 00235 return static_cast<Tuxyz_sgxyzMatrix &> (*(MatrixContainer[MatrixID].MatrixPtr)); 00236 }; 00237 00243 TComplexMatrix & GetComplexMatrix (const TMatrixID MatrixID){ 00244 return static_cast<TComplexMatrix &> (*(MatrixContainer[MatrixID].MatrixPtr)); 00245 }; 00246 00252 TFFTWComplexMatrix & GetFFTWComplexMatrix(const TMatrixID MatrixID){ 00253 return static_cast<TFFTWComplexMatrix &> (*(MatrixContainer[MatrixID].MatrixPtr)); 00254 }; 00255 00261 TLongMatrix & GetLongMatrix (const TMatrixID MatrixID){ 00262 return static_cast<TLongMatrix &> (*(MatrixContainer[MatrixID].MatrixPtr)); 00263 }; 00264 00265 00266 protected: 00267 00268 00269 private: 00270 00272 typedef map<TMatrixID, TMatrixRecord> TMatrixRecordContainer; 00273 00274 00276 TMatrixRecordContainer MatrixContainer; 00277 00279 TMatrixContainer(const TMatrixContainer& src); 00280 00282 TMatrixContainer & operator = (const TMatrixContainer& src); 00283 00285 void PrintErrorAndThrowException(const char * FMT, const string HDF5MatrixName, 00286 const char * File, const int Line); 00287 00288 00289 };// end of TMatrixContainer 00290 00291 00292 #endif /* MATRIXCONTAINER_H */ 00293