kspaceFirstOrder3D-CUDA  1.1
The CUDA/C++ implementation of the k-wave toolbox for the time-domain simulation of acoustic wave fields in 3D
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
KSpaceFirstOrder3DSolver.h
Go to the documentation of this file.
1 /**
2  * @file KSpaceFirstOrder3DSolver.h
3  *
4  * @author Jiri Jaros \n
5  * Faculty of Information Technology \n
6  * Brno University of Technology \n
7  * jarosjir@fit.vutbr.cz
8  *
9  * @brief The header file containing the main class of the project responsible for
10  * the entire 3D fluid simulation.
11  *
12  * @version kspaceFirstOrder3D 3.4
13  *
14  * @date 12 July 2012, 10:27 (created)\n
15  * 10 August 2016, 15:42 (revised)
16  *
17  * @section License
18  * This file is part of the C++ extension of the k-Wave Toolbox
19  * (http://www.k-wave.org).\n Copyright (C) 2016 Jiri Jaros and Bradley Treeby.
20  *
21  * This file is part of the k-Wave. k-Wave is free software: you can redistribute it and/or modify
22  * it under the terms of the GNU Lesser General Public License as published by the Free Software
23  * Foundation, either version 3 of the License, or (at your option) any later version.
24  *
25  * k-Wave is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even
26  * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser
27  * General Public License for more details.
28  *
29  * You should have received a copy of the GNU Lesser General Public License along with k-Wave.
30  * If not, see http://www.gnu.org/licenses/.
31  */
32 
33 #ifndef TKSPACE_FIRST_ORDER_3D_SOLVER_H
34 #define TKSPACE_FIRST_ORDER_3D_SOLVER_H
35 
36 
37 #include <Parameters/Parameters.h>
41 
44 
45 #include <Utils/TimeMeasure.h>
46 
48 
49 
50 /**
51  * @class TKSpaceFirstOrder3DSolver
52  * @brief Class responsible for running the k-space first order 3D method.
53  * @details Class responsible for running the k-space first order 3D method. This class maintain
54  * the whole k-wave (implements the time loop).
55  *
56  */
58 {
59  public:
60  /// Constructor.
62  /// Destructor.
64 
65  /// Memory allocation.
66  virtual void AllocateMemory();
67  /// Memory deallocation.
68  virtual void FreeMemory();
69 
70  /// Load simulation data from the input file.
71  virtual void LoadInputData();
72 
73  /// Compute the k-space simulation.
74  virtual void Compute();
75 
76  /// Get memory usage in MB on the CPU side.
77  size_t GetHostMemoryUsageInMB();
78  /// Get memory usage in MB on the GPU side.
79  size_t GetDeviceMemoryUsageInMB();
80 
81  /// Get code name - release code version.
82  const std::string GetCodeName() const;
83 
84  /// Print the code name and license.
85  void PrintFullNameCodeAndLicense() const;
86 
87  /// Get total simulation time.
88  double GetTotalTime() const;
89  /// Get pre-processing time.
90  double GetPreProcessingTime() const;
91  /// Get data load time.
92  double GetDataLoadTime() const;
93  /// Get simulation time (time loop).
94  double GetSimulationTime() const;
95  /// Get post-processing time.
96  double GetPostProcessingTime() const;
97 
98  /// Get total simulation time cumulated over all legs.
99  double GetCumulatedTotalTime() const;
100  /// Get pre-processing time cumulated over all legs.
101  double GetCumulatedPreProcessingTime() const;
102  /// Get data load time cumulated over all legs.
103  double GetCumulatedDataLoadTime() const;
104  /// Get simulation time (time loop) cumulated over all legs.
105  double GetCumulatedSimulationTime() const;
106  /// Get post-processing time cumulated over all legs.
107  double GetCumulatedPostProcessingTime() const;
108 
109 protected:
110 
111  /// Copy constructor not allowed for public.
113  /// operator = not allowed for public.
115 
116  /// Initialize FFT plans.
117  void InitializeFFTPlans();
118  /// Compute pre-processing phase.
119  void PreProcessingPhase();
120  /// Compute the main time loop of the kspaceFirstOrder3D.
121  void ComputeMainLoop();
122  /// Post processing, and closing the output streams.
123  void PostProcessing();
124 
125  /// Store sensor data.
126  void StoreSensorData();
127  /// Write statistics and header into the output file.
128  void WriteOutputDataInfo();
129  /// Save checkpoint data.
130  void SaveCheckpointData();
131 
132  /// Compute new values of acoustic velocity.
133  void ComputeVelocity();
134  /// Compute new values of acoustic velocity gradients.
136 
137  /// Compute new values of acoustic density for non-linear case.
138  void ComputeDensityNonliner();
139  /// Compute new values of acoustic density for linear case.
140  void ComputeDensityLinear();
141 
142  /// Compute acoustic pressure for nonlinear case
144  /// Compute acoustic pressure for linear case
145  void ComputePressureLinear();
146 
147  /// Add in velocity source
148  void AddVelocitySource();
149  /// Add in pressure source.
150  void AddPressureSource();
151  /// Calculate p0_ ource.
152  void Calculate_p0_source();
153 
154  /// Generate kappa matrix for non-absorbing media.
155  void GenerateKappa();
156  /// Generate kappa matrix, absorb_nabla1, absorb_nabla2 for absorbing media
157  void GenerateKappaAndNablas();
158  /// Generate absorb_tau, absorb_eta for heterogenous media.
159  void GenerateTauAndEta();
160  /// Calculate dt ./ rho0 for non-uniform grids.
161  void GenerateInitialDenisty();
162  /// Calculate square of velocity
163  void Compute_c2();
164 
165  /// Calculate three temporary sums in the new pressure formula, non-linear absorbing case.
167  TRealMatrix& BonA_part,
168  TRealMatrix& du_part);
169 
170  /// Calculate two temporary sums in the new pressure formula, linear absorbing case.
171  void ComputePressurePartsLinear(TRealMatrix& rhoxyz_sum,
172  TRealMatrix& rho0_du_sum);
173 
174  /// Sum sub-terms to calculate new pressure, non-linear case.
175  void SumPressureTermsNonlinear(TRealMatrix& absorb_tau_temp,
176  TRealMatrix& absorb_eta_temp,
177  TRealMatrix& BonA_temp);
178 
179  /// Sum sub-terms to calculate new pressure, linear case.
180  void SumPressureTermsLinear(TRealMatrix& absorb_tau_temp,
181  TRealMatrix& absorb_eta_temp,
182  TRealMatrix& rhoxyz_sum);
183 
184  /// Sum sub-terms for new p, linear lossless case.
186 
187  /// Sum sub-terms for new p, linear lossless case.
189 
190  /// Calculate ux_shifted, uy_shifted and uz_shifted.
192 
193  /// Print progress statistics.
194  void PrintStatistics();
195 
196  /// Is time to checkpoint (save actual state on disk).
197  bool IsTimeToCheckpoint();
198 
199  /// Was the loop interrupted to checkpoint?
200  bool IsCheckpointInterruption() const;
201 
202  /// Check the output file has the correct format and version.
203  void CheckOutputFile();
204 
205  /// Check the checkpoint file has the correct format and version.
206  void CheckCheckpointFile();
207 
208  /// Reads the header of the output file and sets the cumulative elapsed time from the first log.
209  void LoadElapsedTimeFromOutputFile(THDF5_File& o1utputFile);
210 
211  //-------------------------------------- Get matrices ----------------------------------------//
212  /// Get the kappa matrix from the container.
214  {
215  return matrixContainer.GetMatrix<TRealMatrix>(kappa);
216  };
217  /// Get the c^2 matrix from the container.
219  {
221  };
222 
223  /// Get the p matrix from the container.
225  {
227  };
228 
229  /// Get the ux_sgx matrix from the container.
231  {
232  return matrixContainer.GetMatrix<TRealMatrix>(ux_sgx);
233  };
234  /// Get the uy_sgy matrix from the container.
236  {
237  return matrixContainer.GetMatrix<TRealMatrix>(uy_sgy);
238  };
239  /// Get the uz_sgz matrix from the container.
241  {
242  return matrixContainer.GetMatrix<TRealMatrix>(uz_sgz);
243  };
244 
245  /// Get the ux_shifted matrix from the container.
247  {
248  return matrixContainer.GetMatrix<TRealMatrix>(ux_shifted);
249  };
250  /// Get the uy_shifted matrix from the container.
252  {
253  return matrixContainer.GetMatrix<TRealMatrix>(uy_shifted);
254  };
255  /// Get the uz_shifted matrix from the container.
257  {
258  return matrixContainer.GetMatrix<TRealMatrix>(uz_shifted);
259  };
260 
261  /// Get the duxdx matrix from the container.
263  {
264  return matrixContainer.GetMatrix<TRealMatrix>(duxdx);
265  };
266  /// Get the duydy matrix from the container.
268  {
269  return matrixContainer.GetMatrix<TRealMatrix>(duydy);
270  };
271  /// Get the duzdz matrix from the container.
273  {
274  return matrixContainer.GetMatrix<TRealMatrix>(duzdz);
275  };
276 
277  /// Get the dt.*rho0_sgx matrix from the container.
279  {
280  return matrixContainer.GetMatrix<TRealMatrix>(dt_rho0_sgx);
281  };
282  /// Get the dt.*rho0_sgy matrix from the container.
284  {
285  return matrixContainer.GetMatrix<TRealMatrix>(dt_rho0_sgy);
286  };
287  /// Get the dt.*rho0_sgz matrix from the container.
289  {
290  return matrixContainer.GetMatrix<TRealMatrix>(dt_rho0_sgz);
291  };
292 
293  /// Get the rhox matrix from the container.
295  {
296  return matrixContainer.GetMatrix<TRealMatrix>(rhox);
297  };
298  /// Get the rhoy matrix from the container.
300  {
301  return matrixContainer.GetMatrix<TRealMatrix>(rhoy);
302  };
303  /// Get the rhoz matrix from the container.
305  {
306  return matrixContainer.GetMatrix<TRealMatrix>(rhoz);
307  };
308  /// Get the rho0 matrix from the container.
310  {
311  return matrixContainer.GetMatrix<TRealMatrix>(rho0);
312  };
313 
314  /// Get the ddx_k_shift_pos matrix from the container.
316  {
317  return matrixContainer.GetMatrix<TComplexMatrix>(ddx_k_shift_pos);
318  };
319  /// Get the ddy_k_shift_pos matrix from the container.
321  {
322  return matrixContainer.GetMatrix<TComplexMatrix>(ddy_k_shift_pos);
323  };
324  /// Get the ddz_k_shift_pos matrix from the container.
326  {
327  return matrixContainer.GetMatrix<TComplexMatrix>(ddz_k_shift_pos);
328  };
329  /// Get the ddx_k_shift_neg matrix from the container.
331  {
332  return matrixContainer.GetMatrix<TComplexMatrix>(ddx_k_shift_neg);
333  };
334  /// Get the ddy_k_shift_neg matrix from the container.
336  {
337  return matrixContainer.GetMatrix<TComplexMatrix>(ddy_k_shift_neg);
338  };
339  /// Get the ddz_k_shift_neg matrix from the container.
341  {
342  return matrixContainer.GetMatrix<TComplexMatrix>(ddz_k_shift_neg);
343  };
344 
345  /// Get the x_shift_neg_r matrix from the container.
347  {
348  return matrixContainer.GetMatrix<TComplexMatrix>(x_shift_neg_r);
349  };
350  /// Get the y_shift_neg_r from the container.
352  {
353  return matrixContainer.GetMatrix<TComplexMatrix>(y_shift_neg_r);
354  };
355  /// Get the y_shift_neg_r from the container.
357  {
358  return matrixContainer.GetMatrix<TComplexMatrix>(z_shift_neg_r);
359  };
360 
361  /// Get the pml_x_sgx matrix from the container.
363  {
364  return matrixContainer.GetMatrix<TRealMatrix>(pml_x_sgx);
365  };
366  /// Get the pml_y_sgy matrix from the container.
368  {
369  return matrixContainer.GetMatrix<TRealMatrix>(pml_y_sgy);
370  };
371  /// Get the pml_z_sgz matrix from the container.
373  {
374  return matrixContainer.GetMatrix<TRealMatrix>(pml_z_sgz);
375  };
376 
377  /// Get the pml_x matrix from the container.
379  {
380  return matrixContainer.GetMatrix<TRealMatrix>(pml_x);
381  };
382  /// Get the pml_y matrix from the container.
384  {
385  return matrixContainer.GetMatrix<TRealMatrix>(pml_y);
386  };
387  /// Get the pml_z matrix from the container.
389  {
390  return matrixContainer.GetMatrix<TRealMatrix>(pml_z);
391  };
392 
393 
394  /// Get the dxudxn matrix from the container.
396  {
397  return matrixContainer.GetMatrix<TRealMatrix>(dxudxn);
398  };
399  /// Get the dyudyn matrix from the container.
401  {
402  return matrixContainer.GetMatrix<TRealMatrix>(dyudyn);
403  };
404  /// Get the dzudzn matrix from the container.
406  {
407  return matrixContainer.GetMatrix<TRealMatrix>(dzudzn);
408  };
409 
410  /// Get the dxudxn_sgx matrix from the container.
412  {
413  return matrixContainer.GetMatrix<TRealMatrix>(dxudxn_sgx);
414  };
415  /// Get the dyudyn_sgy matrix from the container.
417  {
418  return matrixContainer.GetMatrix<TRealMatrix>(dyudyn_sgy);
419  };
420  /// Get the dzudzn_sgz matrix from the container.
422  {
423  return matrixContainer.GetMatrix<TRealMatrix>(dzudzn_sgz);
424  };
425 
426 
427  /// Get the BonA matrix from the container.
429  {
430  return matrixContainer.GetMatrix<TRealMatrix>(BonA);
431  };
432  /// Get the absorb_tau matrix from the container.
434  {
435  return matrixContainer.GetMatrix<TRealMatrix>(absorb_tau);
436  };
437  /// Get the absorb_eta matrix from the container.
439  {
440  return matrixContainer.GetMatrix<TRealMatrix>(absorb_eta);
441  };
442 
443  /// Get the absorb_nabla1 matrix from the container.
445  {
446  return matrixContainer.GetMatrix<TRealMatrix>(absorb_nabla1);
447  };
448  /// Get the absorb_nabla2 matrix from the container.
450  {
451  return matrixContainer.GetMatrix<TRealMatrix>(absorb_nabla2);
452  };
453 
454 
455  //-- Index matrices --//
456 
457  /// Get the sensor_mask_index matrix from the container.
459  {
460  return matrixContainer.GetMatrix<TIndexMatrix>(sensor_mask_index);
461  };
462 
463 
464  /// Get the sensor_mask_corners matrix from the container.
466  {
467  return matrixContainer.GetMatrix<TIndexMatrix>(sensor_mask_corners);
468  };
469 
470  /// Get the u_source_index matrix from the container.
472  {
473  return matrixContainer.GetMatrix<TIndexMatrix>(u_source_index);
474  };
475  /// Get the p_source_index matrix from the container.
477  {
478  return matrixContainer.GetMatrix<TIndexMatrix>(p_source_index);
479  };
480  /// Get the delay_mask matrix from the container.
482  {
483  return matrixContainer.GetMatrix<TIndexMatrix>(delay_mask);
484  }
485 
486 
487  //-- sources --//
488 
489  /// Get the transducer_source_input matrix from the container.
491  {
492  return matrixContainer.GetMatrix<TRealMatrix>(transducer_source_input);
493  };
494 
495  /// Get the p_source_input matrix from the container.
497  {
498  return matrixContainer.GetMatrix<TRealMatrix>(p_source_input);
499  };
500 
501  /// Get the p0_source_input from the container.
503  {
504  return matrixContainer.GetMatrix<TRealMatrix>(p0_source_input);
505  };
506 
507 
508  /// Get the ux_source_input matrix from the container.
510  {
511  return matrixContainer.GetMatrix<TRealMatrix>(ux_source_input);
512  };
513  /// Get the uy_source_input matrix from the container.
515  {
516  return matrixContainer.GetMatrix<TRealMatrix>(uy_source_input);
517  };
518  /// Get the uz_source_input matrix from the container.
520  {
521  return matrixContainer.GetMatrix<TRealMatrix>(uz_source_input);
522  };
523 
524 
525  //--Temporary matrices --//
526 
527  /// Get the Temp_1_RS3D matrix from the container.
529  {
530  return matrixContainer.GetMatrix<TRealMatrix>(temp_1_real_3D);
531  };
532  /// Get the Temp_2_RS3D matrix from the container.
534  {
535  return matrixContainer.GetMatrix<TRealMatrix>(temp_2_real_3D);
536  };
537  /// Get the Temp_3_RS3D matrix from the container.
539  {
540  return matrixContainer.GetMatrix<TRealMatrix>(temp_3_real_3D);
541  };
542 
543 
544  /// Get the CUFFT_X_temp from the container.
546  {
547  return matrixContainer.GetMatrix<TCUFFTComplexMatrix>(cufft_x_temp);
548  };
549  /// Get the FFT_Y_temp from the container.
551  {
552  return matrixContainer.GetMatrix<TCUFFTComplexMatrix>(cufft_y_temp);
553  };
554  /// Get the FFT_Z_temp from the container.
556  {
557  return matrixContainer.GetMatrix<TCUFFTComplexMatrix>(cufft_z_temp);
558  };
559 
560  /// Get the FFT_shift_temp the container.
562  {
563  return matrixContainer.GetMatrix<TCUFFTComplexMatrix>(cufft_shift_temp);
564  };
565 
566 private:
567 
568  /// Matrix container with all the matrix classes
570  /// Output stream container
572 
573  /// Global parameters of the simulation
575 
576  /// Percentage of the simulation done
577  size_t actPercent;
578 
579  /// This variable is true when calculating first time step after restore from checkpoint (to allow asynchronous IO)
581 
582  /// Total time of the simulation
584  /// Pre-processing time of the simulation
586  /// Data load time of the simulation
588  /// Simulation time of the simulation
590  /// Post-processing time of the simulation
592  /// Iteration time of the simulation
594 
595 };// end of TKSpaceFirstOrder3DSolver
596 //--------------------------------------------------------------------------------------------------
597 
598 #endif /* TKSPACE_FIRST_ORDER_3D_SOLVER_H */
size_t GetDeviceMemoryUsageInMB()
Get memory usage in MB on the GPU side.
TRealMatrix & Get_absorb_nabla1()
Get the absorb_nabla1 matrix from the container.
TRealMatrix & Get_dt_rho0_sgy()
Get the dt.*rho0_sgy matrix from the container.
void ComputePressurePartsLinear(TRealMatrix &rhoxyz_sum, TRealMatrix &rho0_du_sum)
Calculate two temporary sums in the new pressure formula, linear absorbing case.
void AddPressureSource()
Add in pressure source.
double GetTotalTime() const
Get total simulation time.
TRealMatrix & Get_ux_sgx()
Get the ux_sgx matrix from the container.
TRealMatrix & Get_ux_source_input()
Get the ux_source_input matrix from the container.
TRealMatrix & Get_absorb_tau()
Get the absorb_tau matrix from the container.
TComplexMatrix & Get_x_shift_neg_r()
Get the x_shift_neg_r matrix from the container.
void ComputeDensityLinear()
Compute new values of acoustic density for linear case.
void SumPressureLinearLossless()
Sum sub-terms for new p, linear lossless case.
void PrintStatistics()
Print progress statistics.
bool isTimestepRightAfterRestore
This variable is true when calculating first time step after restore from checkpoint (to allow asynch...
TRealMatrix & Get_duydy()
Get the duydy matrix from the container.
TRealMatrix & Get_pml_x()
Get the pml_x matrix from the container.
bool IsCheckpointInterruption() const
Was the loop interrupted to checkpoint?
TRealMatrix & Get_uz_shifted()
Get the uz_shifted matrix from the container.
bool IsTimeToCheckpoint()
Is time to checkpoint (save actual state on disk).
void InitializeFFTPlans()
Initialize FFT plans.
TRealMatrix & Get_dzudzn_sgz()
Get the dzudzn_sgz matrix from the container.
void ComputeGradientVelocity()
Compute new values of acoustic velocity gradients.
void AddVelocitySource()
Add in velocity source.
TComplexMatrix & Get_ddz_k_shift_neg()
Get the ddz_k_shift_neg matrix from the container.
TCUFFTComplexMatrix & Get_cufft_z_temp()
Get the FFT_Z_temp from the container.
TOutputStreamContainer outputStreamContainer
Output stream container.
TRealMatrix & Get_kappa()
Get the kappa matrix from the container.
virtual void FreeMemory()
Memory deallocation.
TTimeMeasure dataLoadTime
Data load time of the simulation.
TRealMatrix & Get_pml_z()
Get the pml_z matrix from the container.
TComplexMatrix & Get_ddz_k_shift_pos()
Get the ddz_k_shift_pos matrix from the container.
TRealMatrix & Get_dxudxn_sgx()
Get the dxudxn_sgx matrix from the container.
TRealMatrix & Get_pml_z_sgz()
Get the pml_z_sgz matrix from the container.
The header file containing the class for real matrices.
TKSpaceFirstOrder3DSolver & operator=(const TKSpaceFirstOrder3DSolver &src)
operator = not allowed for public.
void GenerateKappa()
Generate kappa matrix for non-absorbing media.
TRealMatrix & Get_c2()
Get the c^2 matrix from the container.
TCUFFTComplexMatrix & Get_cufft_y_temp()
Get the FFT_Y_temp from the container.
TRealMatrix & Get_ux_shifted()
Get the ux_shifted matrix from the container.
double GetPreProcessingTime() const
Get pre-processing time.
TTimeMeasure iterationTime
Iteration time of the simulation.
TRealMatrix & Get_dyudyn_sgy()
Get the dyudyn_sgy matrix from the container.
TCUFFTComplexMatrix & Get_cufft_x_temp()
Get the CUFFT_X_temp from the container.
TRealMatrix & Get_p_source_input()
Get the p_source_input matrix from the container.
TTimeMeasure totalTime
Total time of the simulation.
TRealMatrix & Get_temp_3_real_3D()
Get the Temp_3_RS3D matrix from the container.
T & GetMatrix(const TMatrixIdx matrixIdx)
Get the matrix with a specific type from the container.
void Compute_c2()
Calculate square of velocity.
void PrintFullNameCodeAndLicense() const
Print the code name and license.
TRealMatrix & Get_transducer_source_input()
Get the transducer_source_input matrix from the container.
The header file containing the parameters of the simulation.
void SumPressureNonlinearLossless()
Sum sub-terms for new p, linear lossless case.
TComplexMatrix & Get_z_shift_neg_r()
Get the y_shift_neg_r from the container.
const std::string GetCodeName() const
Get code name - release code version.
void CheckCheckpointFile()
Check the checkpoint file has the correct format and version.
TRealMatrix & Get_p0_source_input()
Get the p0_source_input from the container.
TRealMatrix & Get_p()
Get the p matrix from the container.
TComplexMatrix & Get_ddy_k_shift_pos()
Get the ddy_k_shift_pos matrix from the container.
TRealMatrix & Get_BonA()
Get the BonA matrix from the container.
TRealMatrix & Get_dxudxn()
Get the dxudxn matrix from the container.
The header file for class with time measurement.
TRealMatrix & Get_uz_sgz()
Get the uz_sgz matrix from the container.
TMatrixContainer matrixContainer
Matrix container with all the matrix classes.
TRealMatrix & Get_dzudzn()
Get the dzudzn matrix from the container.
double GetSimulationTime() const
Get simulation time (time loop).
void Calculate_p0_source()
Calculate p0_ ource.
TRealMatrix & Get_absorb_nabla2()
Get the absorb_nabla2 matrix from the container.
void GenerateKappaAndNablas()
Generate kappa matrix, absorb_nabla1, absorb_nabla2 for absorbing media.
void StoreSensorData()
Store sensor data.
TRealMatrix & Get_rho0()
Get the rho0 matrix from the container.
void WriteOutputDataInfo()
Write statistics and header into the output file.
TRealMatrix & Get_uz_source_input()
Get the uz_source_input matrix from the container.
void SumPressureTermsLinear(TRealMatrix &absorb_tau_temp, TRealMatrix &absorb_eta_temp, TRealMatrix &rhoxyz_sum)
Sum sub-terms to calculate new pressure, linear case.
TIndexMatrix & Get_sensor_mask_index()
Get the sensor_mask_index matrix from the container.
Class storing all parameters of the simulation.
Definition: Parameters.h:49
TIndexMatrix & Get_delay_mask()
Get the delay_mask matrix from the container.
TRealMatrix & Get_temp_2_real_3D()
Get the Temp_2_RS3D matrix from the container.
The header file containing the class for 64b integer matrices.
void PostProcessing()
Post processing, and closing the output streams.
TRealMatrix & Get_dt_rho0_sgz()
Get the dt.*rho0_sgz matrix from the container.
TIndexMatrix & Get_p_source_index()
Get the p_source_index matrix from the container.
void ComputeMainLoop()
Compute the main time loop of the kspaceFirstOrder3D.
TRealMatrix & Get_dyudyn()
Get the dyudyn matrix from the container.
TRealMatrix & Get_pml_y_sgy()
Get the pml_y_sgy matrix from the container.
TRealMatrix & Get_pml_y()
Get the pml_y matrix from the container.
TIndexMatrix & Get_u_source_index()
Get the u_source_index matrix from the container.
void ComputePressureLinear()
Compute acoustic pressure for linear case.
double GetCumulatedTotalTime() const
Get total simulation time cumulated over all legs.
TIndexMatrix & Get_sensor_mask_corners()
Get the sensor_mask_corners matrix from the container.
TRealMatrix & Get_duxdx()
Get the duxdx matrix from the container.
TCUFFTComplexMatrix & Get_cufft_shift_temp()
Get the FFT_shift_temp the container.
Class implementing the matrix container.
void SumPressureTermsNonlinear(TRealMatrix &absorb_tau_temp, TRealMatrix &absorb_eta_temp, TRealMatrix &BonA_temp)
Sum sub-terms to calculate new pressure, non-linear case.
virtual void LoadInputData()
Load simulation data from the input file.
The header file with the class for complex matrices.
The class for real matrices.
Definition: RealMatrix.h:45
virtual ~TKSpaceFirstOrder3DSolver()
Destructor.
TRealMatrix & Get_dt_rho0_sgx()
Get the dt.*rho0_sgx matrix from the container.
void ComputePressurePartsNonLinear(TRealMatrix &rho_part, TRealMatrix &BonA_part, TRealMatrix &du_part)
Calculate three temporary sums in the new pressure formula, non-linear absorbing case.
double GetCumulatedPostProcessingTime() const
Get post-processing time cumulated over all legs.
The class for 64b unsigned integers (indices). It is used for sensor_mask_index or sensor_corners_mas...
Definition: IndexMatrix.h:46
TTimeMeasure postProcessingTime
Post-processing time of the simulation.
void CheckOutputFile()
Check the output file has the correct format and version.
void GenerateTauAndEta()
Generate absorb_tau, absorb_eta for heterogenous media.
TRealMatrix & Get_uy_sgy()
Get the uy_sgy matrix from the container.
void LoadElapsedTimeFromOutputFile(THDF5_File &o1utputFile)
Reads the header of the output file and sets the cumulative elapsed time from the first log...
TParameters & parameters
Global parameters of the simulation.
Class measuring elapsed time.
Definition: TimeMeasure.h:55
void ComputePressureNonlinear()
Compute acoustic pressure for nonlinear case.
Name space for all CUDA kernels used in the 3D solver.
TTimeMeasure preProcessingTime
Pre-processing time of the simulation.
The header file containing the matrix container and the related matrix record class.
double GetCumulatedSimulationTime() const
Get simulation time (time loop) cumulated over all legs.
size_t actPercent
Percentage of the simulation done.
TRealMatrix & Get_duzdz()
Get the duzdz matrix from the container.
virtual void Compute()
Compute the k-space simulation.
TRealMatrix & Get_rhoz()
Get the rhoz matrix from the container.
TComplexMatrix & Get_ddx_k_shift_neg()
Get the ddx_k_shift_neg matrix from the container.
double GetCumulatedDataLoadTime() const
Get data load time cumulated over all legs.
double GetPostProcessingTime() const
Get post-processing time.
The header file defining the output stream container.
TTimeMeasure simulationTime
Simulation time of the simulation.
TComplexMatrix & Get_ddy_k_shift_neg()
Get the ddy_k_shift_neg matrix from the container.
virtual void AllocateMemory()
Memory allocation.
void SaveCheckpointData()
Save checkpoint data.
void CalculateShiftedVelocity()
Calculate ux_shifted, uy_shifted and uz_shifted.
double GetCumulatedPreProcessingTime() const
Get pre-processing time cumulated over all legs.
TComplexMatrix & Get_ddx_k_shift_pos()
Get the ddx_k_shift_pos matrix from the container.
TRealMatrix & Get_pml_x_sgx()
Get the pml_x_sgx matrix from the container.
double GetDataLoadTime() const
Get data load time.
TRealMatrix & Get_rhox()
Get the rhox matrix from the container.
The class for complex matrices.
Definition: ComplexMatrix.h:54
void ComputeDensityNonliner()
Compute new values of acoustic density for non-linear case.
Class responsible for running the k-space first order 3D method.
TComplexMatrix & Get_y_shift_neg_r()
Get the y_shift_neg_r from the container.
TRealMatrix & Get_temp_1_real_3D()
Get the Temp_1_RS3D matrix from the container.
Class implementing 3D Real-To-Complex and Complex-To-Real transforms using CUDA FFT interface...
TRealMatrix & Get_rhoy()
Get the rhoy matrix from the container.
A container for output streams.
void GenerateInitialDenisty()
Calculate dt ./ rho0 for non-uniform grids.
TRealMatrix & Get_absorb_eta()
Get the absorb_eta matrix from the container.
Class wrapping the HDF5 routines.
Definition: HDF5_File.h:500
void ComputeVelocity()
Compute new values of acoustic velocity.
TRealMatrix & Get_uy_shifted()
Get the uy_shifted matrix from the container.
size_t GetHostMemoryUsageInMB()
Get memory usage in MB on the CPU side.
TRealMatrix & Get_uy_source_input()
Get the uy_source_input matrix from the container.
void PreProcessingPhase()
Compute pre-processing phase.