kspaceFirstOrder3D-OMP 1.0
The C++ implementation of the k-wave toolbox for the time-domain simulation of acoustic wave fields in 3D
MatrixClasses/MatrixContainer.cpp
Go to the documentation of this file.
00001 
00031 #include <stdexcept>
00032 
00033 #include <MatrixClasses/MatrixContainer.h>
00034 
00035 #include <Parameters/Parameters.h>
00036 #include <Utils/ErrorMessages.h>
00037 
00038 //----------------------------------------------------------------------------//
00039 //----------------------------- CONSTANTS ------------------------------------//
00040 //----------------------------------------------------------------------------//
00041 
00042 
00043 
00044 
00045 //============================================================================//
00046 //                              TMatrixRecord                                 //
00047 //============================================================================//
00048 
00049 //----------------------------------------------------------------------------//
00050 //--------------------------- Public methods ---------------------------------//
00051 //----------------------------------------------------------------------------//
00052 
00053 
00058 TMatrixRecord::TMatrixRecord(const TMatrixRecord& src) :
00059     MatrixPtr(src.MatrixPtr), MatrixDataType(src.MatrixDataType),             
00060     DimensionSizes(src.DimensionSizes), LoadData(src.LoadData),  
00061     HDF5MatrixName(src.HDF5MatrixName)
00062 {
00063     
00064 
00065 }// end of TMatrixRecord
00066 //------------------------------------------------------------------------------
00067 
00068 
00074 TMatrixRecord& TMatrixRecord::operator = (const TMatrixRecord& src){
00075 
00076     if (this != &src){
00077         MatrixPtr = src.MatrixPtr;           
00078         MatrixDataType  = src.MatrixDataType;            
00079         DimensionSizes  = src.DimensionSizes;
00080         LoadData        = src.LoadData;
00081         HDF5MatrixName  = src.HDF5MatrixName; 
00082     }
00083     
00084     return *this;
00085    
00086 }// end of operator = 
00087 //------------------------------------------------------------------------------
00088 
00089 
00090     
00099 void TMatrixRecord::SetAllValues(TBaseMatrix *          MatrixPtr,                                                                   
00100                                  const TMatrixDataType  MatrixDataType,                                                                   
00101                                  const TDimensionSizes  DimensionSizes,
00102                                  const bool             LoadData, 
00103                                  const string           HDF5MatrixName){
00104     
00105     this->MatrixPtr        = MatrixPtr;    
00106     this->MatrixDataType   = MatrixDataType;        
00107     this->DimensionSizes   = DimensionSizes;    
00108     this->LoadData         = LoadData;
00109     this->HDF5MatrixName   = HDF5MatrixName; 
00110     
00111 }// end of SetAllValues
00112 //------------------------------------------------------------------------------
00113     
00114 
00115 //----------------------------------------------------------------------------//
00116 //------------------------- Protected methods --------------------------------//
00117 //----------------------------------------------------------------------------//
00118 
00119 //----------------------------------------------------------------------------//
00120 //-------------------------- Private methods ---------------------------------//
00121 //----------------------------------------------------------------------------//
00122 
00123 
00124 
00125 //============================================================================//
00126 //                              TMatrixContainer                              //
00127 //============================================================================//
00128 
00129 //----------------------------------------------------------------------------//
00130 //--------------------------- Public methods ---------------------------------//
00131 //----------------------------------------------------------------------------//
00132 
00133 
00134 
00138 TMatrixContainer::~TMatrixContainer(){
00139     
00140     MatrixContainer.clear();    
00141     
00142 }// end of ~TMatrixContainer
00143 //------------------------------------------------------------------------------
00144 
00145  
00150 void TMatrixContainer::CreateAllObjects(){
00151     
00152     
00153     for (TMatrixRecordContainer::iterator it = MatrixContainer.begin(); it != MatrixContainer.end(); it++){
00154                 
00155         if (it->second.MatrixPtr != NULL) {
00156             PrintErrorAndThrowException(MatrixContainer_ERR_FMT_ReloactaionError, it->second.HDF5MatrixName,
00157                                __FILE__,__LINE__);                    
00158         }                       
00159         
00160         switch (it->second.MatrixDataType){                        
00161             
00162             case TMatrixRecord::mdtReal:  {
00163                  it->second.MatrixPtr =  new TRealMatrix(it->second.DimensionSizes);                    
00164                  break;
00165             }
00166             
00167             case TMatrixRecord::mdtComplex:  {
00168                  it->second.MatrixPtr =  new TComplexMatrix(it->second.DimensionSizes);                    
00169                  break;
00170             }
00171             
00172             case TMatrixRecord::mdtIndex:  {
00173                  it->second.MatrixPtr =  new TLongMatrix(it->second.DimensionSizes);                    
00174                  break;
00175             }
00176             case TMatrixRecord::mdtFFTW:  {
00177                  it->second.MatrixPtr =  new TFFTWComplexMatrix(it->second.DimensionSizes);                    
00178                  break;
00179             }
00180             
00181             case TMatrixRecord::mdtUxyz:  {
00182                  it->second.MatrixPtr =  new Tuxyz_sgxyzMatrix(it->second.DimensionSizes);                    
00183                  break;
00184             }
00185                        
00186             default:{
00187                 PrintErrorAndThrowException(MatrixContainer_ERR_FMT_RecordUnknownDistributionType, it->second.HDF5MatrixName,
00188                                    __FILE__,__LINE__);
00189                                 
00190                 break;    
00191             }    
00192                 
00193         }// switch
00194         
00195     
00196     }// end for
00197        
00198 }// end of CreateAllObjects
00199 //-----------------------------------------------------------------------------
00200     
00201 
00202    
00207 void TMatrixContainer::LoadMatricesDataFromDisk(THDF5_File & HDF5_File){
00208    
00209     for (TMatrixRecordContainer::iterator it = MatrixContainer.begin(); it != MatrixContainer.end(); it++){
00210       
00211         if (it->second.LoadData) {
00212                     
00213             it->second.MatrixPtr->ReadDataFromHDF5File(HDF5_File, it->second.HDF5MatrixName.c_str());        
00214         
00215         }
00216     }       
00217     
00218 }// end of LoadMatricesDataFromDisk
00219 //------------------------------------------------------------------------------
00220 
00221   
00226 void TMatrixContainer::FreeAllMatrices(){
00227      for (TMatrixRecordContainer::iterator it = MatrixContainer.begin(); it != MatrixContainer.end(); it++){
00228         if (it->second.MatrixPtr){
00229           delete it->second.MatrixPtr;
00230           it->second.MatrixPtr = NULL;        
00231         }  
00232     }       
00233      
00234     
00235 }// end of FreeAllMatrices
00236 //------------------------------------------------------------------------------    
00237     
00238 
00239 
00240 
00245 void TMatrixContainer::AddMatricesIntoContainer(){
00246     
00247     TParameters * Params = TParameters::GetInstance();
00248         
00249     TDimensionSizes FullDim = Params->GetFullDimensionSizes();
00250     TDimensionSizes ReducedDim = Params->GetReducedDimensionSizes();
00251     
00252     //----------------------Allocate all matrices ----------------------------//
00253                                              
00254     MatrixContainer[kappa] .SetAllValues(NULL ,TMatrixRecord::mdtReal  , ReducedDim, false, kappa_r_Name);
00255     if (!Params->Get_c0_scalar_flag()) {
00256         MatrixContainer[c2]    .SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , true,  c0_Name);
00257     }
00258     MatrixContainer[p]     .SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, p_Name);
00259      
00260     MatrixContainer[rhox]  .SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, rhox_Name);
00261     MatrixContainer[rhoy]  .SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, rhoy_Name);
00262     MatrixContainer[rhoz]  .SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, rhoz_Name);    
00263     
00264     MatrixContainer[ux_sgx].SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, ux_sgx_Name);
00265     MatrixContainer[uy_sgy].SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, uy_sgy_Name);
00266     MatrixContainer[uz_sgz].SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, uz_sgz_Name);
00267     
00268     MatrixContainer[duxdx] .SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, duxdx_Name);    
00269     MatrixContainer[duydy] .SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, duydy_Name);
00270     MatrixContainer[duzdz] .SetAllValues(NULL ,TMatrixRecord::mdtReal  , FullDim   , false, duzdz_Name);
00271     
00272     if (!Params->Get_rho0_scalar_flag()) {
00273        MatrixContainer[rho0]       .SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim , true, rho0_Name);
00274        MatrixContainer[dt_rho0_sgx].SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim , true, rho0_sgx_Name);
00275        MatrixContainer[dt_rho0_sgy].SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim , true, rho0_sgy_Name);
00276        MatrixContainer[dt_rho0_sgz].SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim , true, rho0_sgz_Name);
00277     }   
00278     
00279 
00280     MatrixContainer[ddx_k_shift_pos].SetAllValues(NULL,TMatrixRecord::mdtComplex, TDimensionSizes(ReducedDim.X, 1, 1), true, ddx_k_shift_pos_r_Name);
00281     MatrixContainer[ddy_k_shift_pos].SetAllValues(NULL,TMatrixRecord::mdtComplex, TDimensionSizes(1, ReducedDim.Y, 1), true, ddy_k_shift_pos_Name);
00282     MatrixContainer[ddz_k_shift_pos].SetAllValues(NULL,TMatrixRecord::mdtComplex, TDimensionSizes(1, 1, ReducedDim.Z), true, ddz_k_shift_pos_Name);
00283 
00284     MatrixContainer[ddx_k_shift_neg].SetAllValues(NULL,TMatrixRecord::mdtComplex, TDimensionSizes(ReducedDim.X ,1, 1), true, ddx_k_shift_neg_r_Name);
00285     MatrixContainer[ddy_k_shift_neg].SetAllValues(NULL,TMatrixRecord::mdtComplex, TDimensionSizes(1, ReducedDim.Y, 1), true, ddy_k_shift_neg_Name);
00286     MatrixContainer[ddz_k_shift_neg].SetAllValues(NULL,TMatrixRecord::mdtComplex, TDimensionSizes(1, 1, ReducedDim.Z), true, ddz_k_shift_neg_Name);
00287 
00288     
00289     MatrixContainer[pml_x_sgx] .SetAllValues(NULL,TMatrixRecord::mdtReal,         TDimensionSizes(FullDim.X, 1, 1), true, pml_x_sgx_Name);
00290     MatrixContainer[pml_y_sgy] .SetAllValues(NULL,TMatrixRecord::mdtReal,         TDimensionSizes(1, FullDim.Y, 1), true, pml_y_sgy_Name);
00291     MatrixContainer[pml_z_sgz] .SetAllValues(NULL,TMatrixRecord::mdtReal,         TDimensionSizes(1, 1, FullDim.Z), true, pml_z_sgz_Name);
00292 
00293     MatrixContainer[pml_x]     .SetAllValues(NULL,TMatrixRecord::mdtReal,         TDimensionSizes(FullDim.X, 1, 1), true, pml_x_Name);
00294     MatrixContainer[pml_y]     .SetAllValues(NULL,TMatrixRecord::mdtReal,         TDimensionSizes(1, FullDim.Y, 1), true, pml_y_Name);
00295     MatrixContainer[pml_z]     .SetAllValues(NULL,TMatrixRecord::mdtReal,         TDimensionSizes(1, 1, FullDim.Z), true, pml_z_Name);
00296     
00297     if (Params->Get_nonlinear_flag()) {
00298         if (! Params->Get_BonA_scalar_flag()) {
00299             MatrixContainer[BonA]      .SetAllValues   (NULL,TMatrixRecord::mdtReal, FullDim   , true, BonA_Name);                                
00300         }
00301     }
00302 
00303     if (Params->Get_absorbing_flag() != 0) {
00304         if (!((Params->Get_c0_scalar_flag()) && (Params->Get_alpha_coeff_scallar_flag()))) {
00305             MatrixContainer[absorb_tau].SetAllValues   (NULL,TMatrixRecord::mdtReal, FullDim   , false, absorb_tau_Name);                                
00306             MatrixContainer[absorb_eta].SetAllValues   (NULL,TMatrixRecord::mdtReal, FullDim   , false, absorb_eta_Name);                                
00307         }
00308         MatrixContainer[absorb_nabla1].SetAllValues(NULL,TMatrixRecord::mdtReal, ReducedDim, false, absorb_nabla1_r_Name);                                
00309         MatrixContainer[absorb_nabla2].SetAllValues(NULL,TMatrixRecord::mdtReal, ReducedDim, false, absorb_nabla2_r_Name
00310         );                                
00311     }
00312         
00313     //-- index matrices --//                
00314     MatrixContainer[sensor_mask_ind].SetAllValues(NULL,TMatrixRecord::mdtIndex, TDimensionSizes(1 ,1, Params->Get_sensor_mask_index_size()), true, sensor_mask_index_Name);                                
00315    
00316     
00317     
00318     // if p0 source flag 
00319     if (Params->Get_p0_source_flag() == 1){
00320        MatrixContainer[p0_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim, true, p0_source_input_Name);                                
00321     }
00322     
00323     
00324     // us_index    
00325     if ((Params->Get_transducer_source_flag() != 0) || 
00326         (Params->Get_ux_source_flag() != 0)         ||
00327         (Params->Get_uy_source_flag() != 0)         ||
00328         (Params->Get_uz_source_flag() != 0)            ){
00329                 
00330                 MatrixContainer[u_source_index].SetAllValues(NULL,TMatrixRecord::mdtIndex,TDimensionSizes(1 ,1, Params->Get_u_source_index_size()), true, u_source_index_Name);                                           
00331     }
00332                 
00333                     
00334         // -- transducer source flag defined
00335     if (Params->Get_transducer_source_flag() != 0) {
00336             
00337         MatrixContainer[delay_mask]             .SetAllValues(NULL,TMatrixRecord::mdtIndex,TDimensionSizes(1 ,1, Params->Get_u_source_index_size())          , true, delay_mask_Name);                                                   
00338         MatrixContainer[transducer_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal ,TDimensionSizes(1 ,1, Params->Get_transducer_source_input_size()), true, transducer_source_input_Name);                                           
00339      
00340     }
00341 
00342     
00343       //-- p variables --// 
00344     if (Params->Get_p_source_flag() != 0){
00345                 
00346            if (Params->Get_p_source_many() == 0) {    // 1D case               
00347                 
00348                MatrixContainer[p_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1 ,1, Params->Get_p_source_flag()), true, p_source_input_Name);                                           
00349                
00350            } else {                                             // 2D case
00351 
00352                MatrixContainer[p_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1 ,Params->Get_p_source_index_size(), Params->Get_p_source_flag()), true, p_source_input_Name);                                           
00353            }   
00354            
00355     
00356            MatrixContainer[p_source_index].SetAllValues(NULL,TMatrixRecord::mdtIndex, TDimensionSizes(1 ,1, Params->Get_p_source_index_size()), true, p_source_index_Name);                                                      
00357     }
00358 
00359     
00360     
00361     //----------------------------uxyz source flags---------------------------//
00362     if (Params->Get_ux_source_flag() != 0) {
00363         if (Params->Get_u_source_many() == 0) { // 1D
00364             
00365             MatrixContainer[ux_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1 ,1, Params->Get_ux_source_flag()), true, ux_source_input_Name);                                           
00366         }
00367         else {                                        // 2D
00368             MatrixContainer[ux_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1 ,Params->Get_u_source_index_size(),Params->Get_ux_source_flag()), true, ux_source_input_Name);                                           
00369             
00370         }                                 
00371     }// ux_source_input
00372     
00373     
00374     if (Params->Get_uy_source_flag() != 0) {
00375         if (Params->Get_u_source_many() == 0) { // 1D
00376             
00377             MatrixContainer[uy_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1 ,1, Params->Get_uy_source_flag()), true, uy_source_input_Name);                                           
00378         }
00379         else {                                        // 2D
00380             MatrixContainer[uy_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1 ,Params->Get_u_source_index_size(),Params->Get_uy_source_flag()), true, uy_source_input_Name);                                           
00381             
00382         }                                 
00383                         
00384     }// uy_source_input
00385     
00386     if (Params->Get_uz_source_flag() != 0) {
00387         if (Params->Get_u_source_many() == 0) { // 1D
00388             
00389             MatrixContainer[uz_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal,TDimensionSizes(1 ,1, Params->Get_uz_source_flag()), true, uz_source_input_Name); 
00390         }
00391         else {                                        // 2D
00392             MatrixContainer[uz_source_input].SetAllValues(NULL,TMatrixRecord::mdtReal,TDimensionSizes(1 ,Params->Get_u_source_index_size(),Params->Get_uz_source_flag()), true, uz_source_input_Name);
00393             
00394         }                                     
00395         
00396         
00397     }// uz_source_input
00398     
00399     
00400     
00401     
00402     //-- Non linear grid
00403     if (Params->Get_nonuniform_grid_flag()!= 0) {            
00404         MatrixContainer[dxudxn].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(FullDim.X, 1, 1), true, dxudxn_Name);                                            
00405         MatrixContainer[dyudyn].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1, FullDim.Y, 1), true, dyudyn_Name);                                            
00406         MatrixContainer[dzudzn].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1 ,1, FullDim.Z), true, dzudzn_Name);                                            
00407 
00408         
00409         MatrixContainer[dxudxn_sgx].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(FullDim.X, 1, 1), true, dxudxn_sgx_Name);                                            
00410         MatrixContainer[dyudyn_sgy].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1, FullDim.Y, 1), true, dyudyn_sgy_Name);
00411         MatrixContainer[dzudzn_sgz].SetAllValues(NULL,TMatrixRecord::mdtReal, TDimensionSizes(1 ,1, FullDim.Z), true, dzudzn_sgz_Name);
00412         
00413     }
00414 
00415     
00416     //--Temporary matrices --//        
00417     // this matrix used to load alpha_coeff for absorb_tau pre-calculation    
00418     
00419     if ((Params->Get_absorbing_flag() != 0) && (!Params->Get_alpha_coeff_scallar_flag())){                
00420         MatrixContainer[Temp_1_RS3D].SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim, true, alpha_coeff_Name);         
00421     }else{
00422         MatrixContainer[Temp_1_RS3D].SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim, false, "");
00423     }
00424 
00425         
00426     MatrixContainer[Temp_2_RS3D].SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim, false, "");                                
00427     MatrixContainer[Temp_3_RS3D].SetAllValues(NULL,TMatrixRecord::mdtReal, FullDim, false, "");                                
00428     
00429     MatrixContainer[FFT_X_temp].SetAllValues(NULL, TMatrixRecord::mdtFFTW, ReducedDim, false, "");                                
00430     MatrixContainer[FFT_Y_temp].SetAllValues(NULL, TMatrixRecord::mdtFFTW, ReducedDim, false, "");                                
00431     MatrixContainer[FFT_Z_temp].SetAllValues(NULL, TMatrixRecord::mdtFFTW, ReducedDim, false, "");                                
00432     
00433     //--- output buffers --//
00434     
00435     
00436     TDimensionSizes SensorDims(Params->Get_sensor_mask_index_size(), 1, 1);
00437             
00438             
00439     if (Params->IsStore_p_rms()){
00440        MatrixContainer[p_sensor_rms].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, p_rms_Name);                                
00441     }
00442         
00443     if (Params->IsStore_p_max()){
00444        MatrixContainer[p_sensor_max].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, p_max_Name);                                
00445     }
00446                 
00447     
00448     if (Params->IsStore_u_rms()){
00449        MatrixContainer[ux_sensor_rms].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, ux_rms_Name);
00450        MatrixContainer[uy_sensor_rms].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, uy_rms_Name);
00451        MatrixContainer[uz_sensor_rms].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, uz_rms_Name);
00452     }
00453         
00454     if (Params->IsStore_u_max()){
00455        MatrixContainer[ux_sensor_max].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, ux_max_Name);
00456        MatrixContainer[uy_sensor_max].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, uy_max_Name);                                
00457        MatrixContainer[uz_sensor_max].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, uz_max_Name);                                       
00458     }
00459             
00460     if (Params->IsStore_I_avg()){
00461        MatrixContainer[Ix_sensor_avg].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, Ix_avg_Name);                                
00462        MatrixContainer[Iy_sensor_avg].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, Iy_avg_Name);                                
00463        MatrixContainer[Iz_sensor_avg].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, Iz_avg_Name);                                
00464        
00465     }
00466         
00467     if (Params->IsStore_I_max()){
00468        MatrixContainer[Ix_sensor_max].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, Ix_max_Name);
00469        MatrixContainer[Iy_sensor_max].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, Iy_max_Name);
00470        MatrixContainer[Iz_sensor_max].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, Iz_max_Name);
00471     }
00472     
00473     
00474     // in case of intensity create buffer for time staggered values
00475     if ((Params->IsStore_I_avg()) || (Params->IsStore_I_max())){        
00476        MatrixContainer[p_sensor_i_1_raw] .SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, "");
00477        MatrixContainer[ux_sensor_i_1_agr_2].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, "");
00478        MatrixContainer[uy_sensor_i_1_agr_2].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, "");
00479        MatrixContainer[uz_sensor_i_1_agr_2].SetAllValues(NULL,TMatrixRecord::mdtReal, SensorDims, false, "");        
00480     }
00481     
00482     
00483     
00484 }// end of InitMatrixContainer
00485 //------------------------------------------------------------------------------
00486 
00487 
00488 //----------------------------------------------------------------------------//
00489 //------------------------- Protected methods --------------------------------//
00490 //----------------------------------------------------------------------------//
00491 
00492 
00493 //----------------------------------------------------------------------------//
00494 //-------------------------- Private methods ---------------------------------//
00495 //----------------------------------------------------------------------------//
00496     
00497 
00508 void TMatrixContainer::PrintErrorAndThrowException(const char * FMT, const string HDF5MatrixName,           
00509                                           const char * File, const int Line){
00510     
00511      fprintf(stderr,FMT, HDF5MatrixName.c_str(), File, Line);                                
00512      throw bad_alloc();
00513     
00514 }// end of PrintErrorAndAbort
00515 //------------------------------------------------------------------------------    
00516     
00517 
00518 
 All Classes Files Functions Variables Typedefs Enumerations