kspaceFirstOrder3D-OMP  1.1
The C++ implementation of the k-wave toolbox for the time-domain simulation of acoustic wave fields in 3D
 All Classes Files Functions Variables Typedefs Enumerations Friends Pages
KSpaceFirstOrder3DSolver.h
Go to the documentation of this file.
1 /**
2  * @file KSpaceFirstOrder3DSolver.h
3  * @author Jiri Jaros \n
4  * Faculty of Information Technology\n
5  * Brno University of Technology \n
6  * jarosjir@fit.vutbr.cz
7  *
8  * @brief The header file containing the main class of the project
9  * responsible for the entire simulation.
10  *
11  * @version kspaceFirstOrder3D 2.15
12  *
13  * @date 12 July 2012, 10:27 (created)\n
14  * 29 September 2014, 14:45 (revised)
15  *
16  * @section License
17  * This file is part of the C++ extension of the k-Wave Toolbox (http://www.k-wave.org).\n
18  * Copyright (C) 2014 Jiri Jaros and Bradley Treeby.
19  *
20  * This file is part of k-Wave. k-Wave is free software: you can redistribute it
21  * and/or modify it under the terms of the GNU Lesser General Public License as
22  * published by the Free Software Foundation, either version 3 of the License,
23  * or (at your option) any later version.
24  *
25  * k-Wave is distributed in the hope that it will be useful, but
26  * WITHOUT ANY WARRANTY; without even the implied warranty of
27  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
28  * See the GNU Lesser General Public License for more details.
29  *
30  * You should have received a copy of the GNU Lesser General Public License
31  * along with k-Wave. If not, see <http://www.gnu.org/licenses/>.
32  */
33 
34 #ifndef TKSPACE3DSOLVER_H
35 #define TKSPACE3DSOLVER_H
36 
37 #include <fftw3.h>
38 
39 #include <Parameters/Parameters.h>
40 
48 
49 #include <Utils/TimeMeasure.h>
50 
51 using namespace std;
52 
53 
54 /**
55  * @class TKSpaceFirstOrder3DSolver
56  * @brief Class responsible for running the k-space first order 3D method.
57  * @details Class responsible for running the k-space first order 3D method. This
58  * class maintain the whole k-wave (implements the time loop).
59  *
60  */
62 {
63  public:
64  /// Constructor.
66  /// Destructor.
67  virtual ~TKSpaceFirstOrder3DSolver();
68 
69 
70  /// Memory allocation.
71  virtual void AllocateMemory();
72  /// Memory deallocation.
73  virtual void FreeMemory();
74 
75  /// Load simulation data from the input file.
76  virtual void LoadInputData();
77 
78  /// Compute the 3D kspace first order simulation.
79  virtual void Compute();
80 
81  /// Print parameters of the simulation.
82  virtual void PrintParametersOfSimulation(FILE * file);
83 
84  /// Get memory usage in MB.
85  virtual size_t ShowMemoryUsageInMB();
86 
87  /// Get code name.
88  string GetCodeName() {return "kspaceFirstOrder3D-OMP v1.1"; };
89 
90  /// Print the code name and license.
91  void PrintFullNameCodeAndLicense(FILE * file);
92 
93  /// Set processor affinity.
94  void SetProcessorAffinity();
95 
96  /// Get total simulation time.
97  double GetTotalTime() const { return TotalTime.GetElapsedTime(); };
98 
99  /// Get pre-processing time.
100  double GetPreProcessingTime() const { return PreProcessingTime.GetElapsedTime(); };
101 
102  /// Get data load time.
103  double GetDataLoadTime() const { return DataLoadTime.GetElapsedTime(); };
104 
105  /// Get simulation time (time loop).
106  double GetSimulationTime() const { return SimulationTime.GetElapsedTime(); };
107 
108  /// Get post-processing time.
109  double GetPostProcessingTime() const { return PostProcessingTime.GetElapsedTime();};
110 
111  /// Get total simulation time cumulated over all legs.
112  double GetCumulatedTotalTime() const { return TotalTime.GetCumulatedElapsedTimeOverAllLegs(); };
113  /// Get pre-processing time cumulated over all legs.
114  double GetCumulatedPreProcessingTime() const { return PreProcessingTime.GetCumulatedElapsedTimeOverAllLegs(); };
115  /// Get data load time cumulated over all legs.
116  double GetCumulatedDataLoadTime() const { return DataLoadTime.GetCumulatedElapsedTimeOverAllLegs(); };
117  /// Get simulation time (time loop) cumulated over all legs.
118  double GetCumulatedSimulationTime() const { return SimulationTime.GetCumulatedElapsedTimeOverAllLegs(); };
119  /// Get post-processing time cumulated over all legs.
120  double GetCumulatedPostProcessingTime() const { return PostProcessingTime.GetCumulatedElapsedTimeOverAllLegs(); };
121 
122  protected:
123 
124  /// Copy constructor not allowed for public.
126  /// operator = not allowed for public.
127  TKSpaceFirstOrder3DSolver& operator = (const TKSpaceFirstOrder3DSolver& src);
128 
129  /// Initialize FFT plans.
130  void InitializeFFTWPlans();
131 
132  /// Compute pre-processing phase.
133  void PreProcessingPhase();
134 
135  /// Compute the main time loop of the kspaceFirstOrder3D.
136  void Compute_MainLoop();
137 
138  /// Post processing, and closing the output streams.
139  void PostPorcessing();
140 
141  /// Store sensor data.
142  void StoreSensorData();
143 
144  /// Save checkpoint data.
145  void SaveCheckpointData();
146 
147  /// Write statistics and header into the output file.
148  void WriteOutputDataInfo();
149 
150 
151  /// compute new values of for ux_sgx, uy_sgy, uz_sgz.
152  void Compute_uxyz();
153 
154  /// Compute new values of for duxdx, duydy, dzdz.
155  void Compute_duxyz();
156 
157 
158  /// Compute new values of rhox, rhoy, rhoz for non-linear case.
159  void Compute_rhoxyz_nonlinear();
160  /// Compute new values of rhox, rhoy, rhoz for linear case.
161  void Compute_rhoxyz_linear();
162 
163  /// Add u source to the particle velocity.
164  void Add_u_source();
165 
166  /// Add in pressure source.
167  void Add_p_source();
168 
169 
170  /// Generate kappa matrix for non-absorbing media.
171  void Generate_kappa();
172 
173  /// Generate kappa matrix, absorb_nabla1, absorb_nabla2 for absorbing media.
174  void Generate_kappa_absorb_nabla1_absorb_nabla2();
175 
176  /// Generate absorb_tau, absorb_eta for heterogenous media.
177  void Generate_absorb_tau_absorb_eta_matrix();
178 
179  /// Calculate dt ./ rho0 for non-uniform grids.
180  void Calculate_dt_rho0_non_uniform();
181 
182  /// Calculate p0_source.
183  void Calculate_p0_source();
184 
185  /// Calculate c^2.
186  void Compute_c2();
187 
188  /// Compute part of the new velocity - gradient of p.
189  void Compute_ddx_kappa_fft_p(TRealMatrix& X_Matrix,
190  TFFTWComplexMatrix& FFT_X,
191  TFFTWComplexMatrix& FFT_Y,
192  TFFTWComplexMatrix& FFT_Z,
193  const TRealMatrix& kappa,
194  const TComplexMatrix& ddx,
195  const TComplexMatrix& ddy,
196  const TComplexMatrix& ddz);
197 
198  /// Calculate new p, non-linear case.
199  void Compute_new_p_nonlinear();
200  /// Calculate new p linear-case, absorbing.
201  void Compute_new_p_linear();
202 
203 
204  /// Calculate three temporary sums in the new pressure formula, non-linear absorbing case, SSE2 version.
205  void Calculate_SumRho_BonA_SumDu_SSE2(TRealMatrix& RHO_Temp,
206  TRealMatrix& BonA_Temp,
207  TRealMatrix& Sum_du);
208  /// Calculate two temporary sums in the new pressure formula, linear absorbing case.
209  void Calculate_SumRho_SumRhoDu(TRealMatrix& Sum_rhoxyz,
210  TRealMatrix& Sum_rho0_du);
211 
212  /// Compute absorbing term with abosrb_nabla1 and absorb_nabla2, SSE2 version.
213  void Compute_Absorb_nabla1_2_SSE2(TFFTWComplexMatrix& FFT_1,
214  TFFTWComplexMatrix& FFT_2);
215 
216  /// Sum sub-terms to calculate new pressure, non-linear case.
217  void Sum_Subterms_nonlinear(TRealMatrix& Absorb_tau_temp,
218  TRealMatrix& Absorb_eta_temp,
219  TRealMatrix& BonA_temp);
220 
221  /// Sum sub-terms to calculate new pressure, linear case.
222  void Sum_Subterms_linear(TRealMatrix& Absorb_tau_temp,
223  TRealMatrix& Absorb_eta_temp,
224  TRealMatrix& Sum_rhoxyz);
225 
226  /// Sum sub-terms for new p, linear lossless case.
227  void Sum_new_p_nonlinear_lossless();
228  /// Sum sub-terms for new p, linear lossless case.
229  void Sum_new_p_linear_lossless();
230 
231  /// Calculate ux_shifted, uy_shifted and uz_shifted.
232  void Calculate_shifted_velocity();
233 
234  /// Print progress statistics.
235  void PrintStatisitcs();
236 
237  /// Print the header of the progress statistics.
238  void PrintOtputHeader();
239 
240  /// Is time to checkpoint (save actual state on disk)?
241  bool IsTimeToCheckpoint();
242 
243  /// Was the loop interrupted to checkpoint?
245  {
246  return (Parameters->Get_t_index() != Parameters->Get_Nt());
247  };
248 
249  /// Check the output file has the correct format and version.
250  void CheckOutputFile();
251  /// Reads the header of the output file and sets the cumulative elapsed time from the first log.
252  void RestoreCumulatedElapsedFromOutputFile();
253 
254  /// Check the checkpoint file has the correct format and version.
255  void CheckCheckpointFile();
256 
257 
258  //------------------------- Get matrices --------------------------------//
259  /// Get the kappa matrix from the container.
261  {
262  return MatrixContainer.GetMatrix<TRealMatrix>(kappa);
263  };
264  /// Get the c^2 matrix from the container.
266  {
267  return MatrixContainer.GetMatrix<TRealMatrix>(c2);
268  };
269 
270  /// Get the p matrix from the container.
272  {
273  return MatrixContainer.GetMatrix<TRealMatrix>(p);
274  };
275 
276  /// Get the ux_sgx matrix from the container.
278  {
279  return MatrixContainer.GetMatrix<Tuxyz_sgxyzMatrix>(ux_sgx);
280  };
281  /// Get the uy_sgy matrix from the container.
283  {
284  return MatrixContainer.GetMatrix<Tuxyz_sgxyzMatrix>(uy_sgy);
285  };
286  /// Get the uz_sgz matrix from the container.
288  {
289  return MatrixContainer.GetMatrix<Tuxyz_sgxyzMatrix>(uz_sgz);
290  };
291 
292  /// Get the ux_shifted matrix from the container.
294  {
295  return MatrixContainer.GetMatrix<Tuxyz_sgxyzMatrix>(ux_shifted);
296  };
297  /// Get the uy_shifted matrix from the container.
299  {
300  return MatrixContainer.GetMatrix<Tuxyz_sgxyzMatrix>(uy_shifted);
301  };
302  /// Get the uz_shifted matrix from the container.
304  {
305  return MatrixContainer.GetMatrix<Tuxyz_sgxyzMatrix>(uz_shifted);
306  };
307 
308  /// Get the duxdx matrix from the container.
310  {
311  return MatrixContainer.GetMatrix<TRealMatrix>(duxdx);
312  };
313  /// Get the duydy matrix from the container.
315  {
316  return MatrixContainer.GetMatrix<TRealMatrix>(duydy);
317  };
318  /// Get the duzdz matrix from the container.
320  {
321  return MatrixContainer.GetMatrix<TRealMatrix>(duzdz);
322  };
323 
324  /// Get the dt.*rho0_sgx matrix from the container.
326  {
327  return MatrixContainer.GetMatrix<TRealMatrix>(dt_rho0_sgx);
328  };
329  /// Get the dt.*rho0_sgy matrix from the container.
331  {
332  return MatrixContainer.GetMatrix<TRealMatrix>(dt_rho0_sgy);
333  };
334  /// Get the dt.*rho0_sgz matrix from the container.
336  {
337  return MatrixContainer.GetMatrix<TRealMatrix>(dt_rho0_sgz);
338  };
339 
340  /// Get the rhox matrix from the container.
342  {
343  return MatrixContainer.GetMatrix<TRealMatrix>(rhox);
344  };
345  /// Get the rhoy matrix from the container.
347  {
348  return MatrixContainer.GetMatrix<TRealMatrix>(rhoy);
349  };
350  /// Get the rhoz matrix from the container.
352  {
353  return MatrixContainer.GetMatrix<TRealMatrix>(rhoz);
354  };
355  /// Get the rho0 matrix from the container.
357  {
358  return MatrixContainer.GetMatrix<TRealMatrix>(rho0);
359  };
360 
361  /// Get the ddx_k_shift_pos matrix from the container.
363  {
364  return MatrixContainer.GetMatrix<TComplexMatrix>(ddx_k_shift_pos);
365  };
366  /// Get the ddy_k_shift_pos matrix from the container.
368  {
369  return MatrixContainer.GetMatrix<TComplexMatrix>(ddy_k_shift_pos);
370  };
371  /// Get the ddz_k_shift_pos matrix from the container.
373  {
374  return MatrixContainer.GetMatrix<TComplexMatrix>(ddz_k_shift_pos);
375  };
376  /// Get the ddx_k_shift_neg matrix from the container.
378  {
379  return MatrixContainer.GetMatrix<TComplexMatrix>(ddx_k_shift_neg);
380  };
381  /// Get the ddy_k_shift_neg matrix from the container.
383  {
384  return MatrixContainer.GetMatrix<TComplexMatrix>(ddy_k_shift_neg);
385  };
386  /// Get the ddz_k_shift_neg matrix from the container.
388  {
389  return MatrixContainer.GetMatrix<TComplexMatrix>(ddz_k_shift_neg);
390  };
391 
392  /// Get the x_shift_neg_r matrix from the container.
394  {
395  return MatrixContainer.GetMatrix<TComplexMatrix>(x_shift_neg_r);
396  };
397  /// Get the y_shift_neg_r from the container.
399  {
400  return MatrixContainer.GetMatrix<TComplexMatrix>(y_shift_neg_r);
401  };
402  /// Get the y_shift_neg_r from the container.
404  {
405  return MatrixContainer.GetMatrix<TComplexMatrix>(z_shift_neg_r);
406  };
407 
408  /// Get the pml_x_sgx matrix from the container.
410  {
411  return MatrixContainer.GetMatrix<TRealMatrix>(pml_x_sgx);
412  };
413  /// Get the pml_y_sgy matrix from the container.
415  {
416  return MatrixContainer.GetMatrix<TRealMatrix>(pml_y_sgy);
417  };
418  /// Get the pml_z_sgz matrix from the container.
420  {
421  return MatrixContainer.GetMatrix<TRealMatrix>(pml_z_sgz);
422  };
423 
424  /// Get the pml_x matrix from the container.
426  {
427  return MatrixContainer.GetMatrix<TRealMatrix>(pml_x);
428  };
429  /// Get the pml_y matrix from the container.
431  {
432  return MatrixContainer.GetMatrix<TRealMatrix>(pml_y);
433  };
434  /// Get the pml_z matrix from the container.
436  {
437  return MatrixContainer.GetMatrix<TRealMatrix>(pml_z);
438  };
439 
440 
441  /// Get the dxudxn matrix from the container.
443  {
444  return MatrixContainer.GetMatrix<TRealMatrix>(dxudxn);
445  };
446  /// Get the dyudyn matrix from the container.
448  {
449  return MatrixContainer.GetMatrix<TRealMatrix>(dyudyn);
450  };
451  /// Get the dzudzn matrix from the container.
453  {
454  return MatrixContainer.GetMatrix<TRealMatrix>(dzudzn);
455  };
456 
457  /// Get the dxudxn_sgx matrix from the container.
459  {
460  return MatrixContainer.GetMatrix<TRealMatrix>(dxudxn_sgx);
461  };
462  /// Get the dyudyn_sgy matrix from the container.
464  {
465  return MatrixContainer.GetMatrix<TRealMatrix>(dyudyn_sgy);
466  };
467  /// Get the dzudzn_sgz matrix from the container.
469  {
470  return MatrixContainer.GetMatrix<TRealMatrix>(dzudzn_sgz);
471  };
472 
473 
474  /// Get the BonA matrix from the container.
476  {
477  return MatrixContainer.GetMatrix<TRealMatrix>(BonA);
478  };
479  /// Get the absorb_tau matrix from the container.
481  {
482  return MatrixContainer.GetMatrix<TRealMatrix>(absorb_tau);
483  };
484  /// Get the absorb_eta matrix from the container.
486  {
487  return MatrixContainer.GetMatrix<TRealMatrix>(absorb_eta);
488  };
489 
490  /// Get the absorb_nabla1 matrix from the container.
492  {
493  return MatrixContainer.GetMatrix<TRealMatrix>(absorb_nabla1);
494  };
495  /// Get the absorb_nabla2 matrix from the container.
497  {
498  return MatrixContainer.GetMatrix<TRealMatrix>(absorb_nabla2);
499  };
500 
501 
502  //-- Index matrices --//
503 
504  /// Get the sensor_mask_index matrix from the container.
506  {
507  return MatrixContainer.GetMatrix<TIndexMatrix>(sensor_mask_index);
508  };
509  /// Get the sensor_mask_corners matrix from the container.
511  {
512  return MatrixContainer.GetMatrix<TIndexMatrix>(sensor_mask_corners);
513  };
514 
515  /// Get the u_source_index matrix from the container.
517  {
518  return MatrixContainer.GetMatrix<TIndexMatrix>(u_source_index);
519  };
520  /// Get the p_source_index matrix from the container.
522  {
523  return MatrixContainer.GetMatrix<TIndexMatrix>(p_source_index);
524  };
525  /// Get the delay_mask matrix from the container.
527  {
528  return MatrixContainer.GetMatrix<TIndexMatrix>(delay_mask);
529  }
530 
531 
532  //-- sources --//
533 
534  /// Get the transducer_source_input matrix from the container.
536  {
537  return MatrixContainer.GetMatrix<TRealMatrix>(transducer_source_input);
538  };
539 
540  /// Get the p_source_input matrix from the container.
542  {
543  return MatrixContainer.GetMatrix<TRealMatrix>(p_source_input);
544  };
545 
546  /// Get the p0_source_input from the container.
548  {
549  return MatrixContainer.GetMatrix<TRealMatrix>(p0_source_input);
550  };
551 
552 
553  /// Get the ux_source_input matrix from the container.
555  {
556  return MatrixContainer.GetMatrix<TRealMatrix>(ux_source_input);
557  };
558  /// Get the uy_source_input matrix from the container.
560  {
561  return MatrixContainer.GetMatrix<TRealMatrix>(uy_source_input);
562  };
563  /// Get the uz_source_input matrix from the container.
565  {
566  return MatrixContainer.GetMatrix<TRealMatrix>(uz_source_input);
567  };
568 
569 
570  //--Temporary matrices --//
571 
572  /// Get the Temp_1_RS3D matrix from the container.
574  {
575  return MatrixContainer.GetMatrix<TRealMatrix>(Temp_1_RS3D);
576  };
577  /// Get the Temp_2_RS3D matrix from the container.
579  {
580  return MatrixContainer.GetMatrix<TRealMatrix>(Temp_2_RS3D);
581  };
582  /// Get the Temp_3_RS3D matrix from the container.
584  {
585  return MatrixContainer.GetMatrix<TRealMatrix>(Temp_3_RS3D);
586  };
587 
588 
589  /// Get the FFT_X_temp from the container.
591  {
592  return MatrixContainer.GetMatrix<TFFTWComplexMatrix>(FFT_X_temp);
593  };
594  /// Get the FFT_Y_temp from the container.
596  {
597  return MatrixContainer.GetMatrix<TFFTWComplexMatrix>(FFT_Y_temp);
598  };
599  /// Get the FFT_Z_temp from the container.
601  {
602  return MatrixContainer.GetMatrix<TFFTWComplexMatrix>(FFT_Z_temp);
603  };
604  /// Get the FFT_shift_temp the container.
606  {
607  return MatrixContainer.GetMatrix<TFFTWComplexMatrix>(FFT_shift_temp);
608  };
609 
610  private:
611 
612  /// Matrix container with all the matrix classes.
613  TMatrixContainer MatrixContainer;
614  /// Output stream container.
616 
617  /// Percentage of the simulation done.
618  size_t ActPercent;
619 
620  /// Global parameters of the simulation.
622 
623  /// Total time of the simulation.
625  /// Pre-processing time of the simulation.
627  /// Data load time of the simulation.
629  /// Simulation time of the simulation.
631  /// Post-processing time of the simulation.
633  /// Iteration time of the simulation.
635 
636 
637 };// end of TKSpace3DSolver
638 //------------------------------------------------------------------------------
639 
640 #endif /* TKSPACE3DSOLVER_H */
641 
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.
TTimeMeasure TotalTime
Total time of the simulation.
double GetTotalTime() const
Get total simulation time.
TFFTWComplexMatrix & Get_FFT_X_temp()
Get the FFT_X_temp from the container.
TRealMatrix & Get_ux_source_input()
Get the ux_source_input matrix from the container.
TRealMatrix & Get_Temp_3_RS3D()
Get the Temp_3_RS3D 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.
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_dzudzn_sgz()
Get the dzudzn_sgz matrix from the container.
TComplexMatrix & Get_ddz_k_shift_neg()
Get the ddz_k_shift_neg matrix from the container.
TRealMatrix & Get_kappa()
Get the kappa matrix from the container.
TTimeMeasure SimulationTime
Simulation time of the simulation.
TTimeMeasure PostProcessingTime
Post-processing 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.
TOutputStreamContainer OutputStreamContainer
Output stream container.
Tuxyz_sgxyzMatrix & Get_ux_sgx()
Get the ux_sgx matrix from the container.
TParameters * Parameters
Global parameters of the simulation.
TRealMatrix & Get_c2()
Get the c^2 matrix from the container.
Tuxyz_sgxyzMatrix & Get_uy_sgy()
Get the uy_sgy matrix from the container.
The header file containing the particle velocity matrix.
Tuxyz_sgxyzMatrix & Get_uy_shifted()
Get the uy_shifted matrix from the container.
double GetPreProcessingTime() const
Get pre-processing time.
TRealMatrix & Get_dyudyn_sgy()
Get the dyudyn_sgy matrix from the container.
TRealMatrix & Get_p_source_input()
Get the p_source_input matrix from the container.
TRealMatrix & Get_transducer_source_input()
Get the transducer_source_input matrix from the container.
The header file containing the parameters of the simulation.
TComplexMatrix & Get_z_shift_neg_r()
Get the y_shift_neg_r from the container.
Tuxyz_sgxyzMatrix & Get_uz_sgz()
Get the uz_sgz matrix from the container.
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.
TTimeMeasure DataLoadTime
Data load time of the simulation.
TRealMatrix & Get_BonA()
Get the BonA matrix from the container.
TRealMatrix & Get_dxudxn()
Get the dxudxn matrix from the container.
The header file containing the class measuring elapsed time.
Class implementing 3D Real-To-Complex and Complex-To-Real transforms using FFTW interface.
TRealMatrix & Get_dzudzn()
Get the dzudzn matrix from the container.
double GetSimulationTime() const
Get simulation time (time loop).
TRealMatrix & Get_absorb_nabla2()
Get the absorb_nabla2 matrix from the container.
TRealMatrix & Get_rho0()
Get the rho0 matrix from the container.
TRealMatrix & Get_uz_source_input()
Get the uz_source_input matrix from the container.
The header file of classes responsible for storing output quantities into the output HDF5 file...
Tuxyz_sgxyzMatrix & Get_uz_shifted()
Get the uz_shifted matrix from the container.
TIndexMatrix & Get_sensor_mask_index()
Get the sensor_mask_index matrix from the container.
TFFTWComplexMatrix & Get_FFT_shift_temp()
Get the FFT_shift_temp the container.
Class storing all parameters of the simulation.
Definition: Parameters.h:51
TIndexMatrix & Get_delay_mask()
Get the delay_mask matrix from the container.
The header file containing the class for 64b integer matrices.
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.
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.
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.
Class implementing the matrix container.
TRealMatrix & Get_Temp_1_RS3D()
Get the Temp_1_RS3D matrix from the container.
The header file with the class for complex matrices.
TFFTWComplexMatrix & Get_FFT_Z_temp()
Get the FFT_Z_temp from the container.
The class for real matrices.
Definition: RealMatrix.h:48
TRealMatrix & Get_dt_rho0_sgx()
Get the dt.*rho0_sgx matrix from the container.
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:50
Class measuring elapsed time.
Definition: TimeMeasure.h:51
TTimeMeasure PreProcessingTime
Pre-processing time of the simulation.
string GetCodeName()
Get code name.
The header file containing the matrix container.
double GetCumulatedSimulationTime() const
Get simulation time (time loop) cumulated over all legs.
TFFTWComplexMatrix & Get_FFT_Y_temp()
Get the FFT_Y_temp from the container.
TTimeMeasure IterationTime
Iteration time of the simulation.
TRealMatrix & Get_duzdz()
Get the duzdz matrix from the container.
TRealMatrix & Get_Temp_2_RS3D()
Get the Temp_2_RS3D matrix from the container.
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.
TComplexMatrix & Get_ddy_k_shift_neg()
Get the ddy_k_shift_neg matrix from the container.
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.
Tuxyz_sgxyzMatrix & Get_ux_shifted()
Get the ux_shifted 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 velocity matrix.
The class for complex matrices.
Definition: ComplexMatrix.h:64
size_t ActPercent
Percentage of the simulation done.
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_rhoy()
Get the rhoy matrix from the container.
A container for output streams.
TRealMatrix & Get_absorb_eta()
Get the absorb_eta matrix from the container.
The header file containing the class that implements 3D FFT using the FFTW interface.
TRealMatrix & Get_uy_source_input()
Get the uy_source_input matrix from the container.