WarpX
BTDiagnostics.H
Go to the documentation of this file.
1 /* Copyright 2021 Revathi Jambunathan
2  *
3  * This file is part of WarpX.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 #ifndef WARPX_BTDIAGNOSTICS_H_
8 #define WARPX_BTDIAGNOSTICS_H_
9 
10 #include "Diagnostics.H"
12 #include "Utils/WarpXConst.H"
14 
15 #include <AMReX_Box.H>
16 #include <AMReX_Geometry.H>
17 #include <AMReX_IntVect.H>
18 #include <AMReX_MultiFab.H>
19 #include <AMReX_REAL.H>
20 #include <AMReX_RealBox.H>
21 #include <AMReX_Vector.H>
22 
23 #include <limits>
24 #include <memory>
25 #include <string>
26 
27 class BTDiagnostics final : public Diagnostics
28 {
29 public:
30 
31  BTDiagnostics (int i, const std::string& name);
32 
33 private:
35  bool m_plot_raw_fields = false;
39  void ReadParameters ();
49  void Flush (int i_buffer, bool force_flush) override;
58  bool DoDump (int step, int i_buffer, bool force_flush=false) override;
66  bool DoComputeAndPack (int step, bool force_flush=false) override;
68  void DerivedInitData () override;
75  void InitializeFieldFunctors (int lev) override;
82  void InitializeFieldFunctorsRZopenPMD (int lev) override;
94  void AddRZModesToOutputNames (const std::string& field, int ncomp, bool cellcenter_data);
99  void InitializeParticleBuffer () override;
103  void PrepareBufferData () override;
105  void UpdateBufferData () override;
113  void PrepareFieldDataForOutput () override;
115  void PrepareParticleDataForOutput() override;
116 
125  bool GetZSliceInDomainFlag (int i_buffer, int lev);
126 
133  bool GetKIndexInSnapshotBoxFlag (int i_buffer, int lev);
134 
140  void InitializeBufferData ( int i_buffer , int lev, bool restart=false) override;
149 
157  amrex::Real m_gamma_boost;
158  amrex::Real m_beta_boost;
164 
166  int m_num_snapshots_lab = std::numeric_limits<int>::lowest();
168  amrex::Real m_dt_snapshots_lab = std::numeric_limits<amrex::Real>::lowest();
172  amrex::Real m_dz_snapshots_lab = 0.0;
173 
175  int m_buffer_size = 256;
176 
254  void DefineCellCenteredMultiFab(int lev);
261  void DefineFieldBufferMultiFab (int i_buffer, int lev);
262 
269  void DefineSnapshotGeometry (int i_buffer, int lev);
270 
275  [[nodiscard]] amrex::Real UpdateCurrentZBoostCoordinate(amrex::Real t_lab, amrex::Real t_boost) const
276  {
277  const amrex::Real current_z_boost = (t_lab / m_gamma_boost - t_boost) * PhysConst::c / m_beta_boost;
278  return current_z_boost;
279  }
284  [[nodiscard]] amrex::Real UpdateCurrentZLabCoordinate(amrex::Real t_lab, amrex::Real t_boost) const
285  {
286  const amrex::Real current_z_lab = (t_lab - t_boost / m_gamma_boost ) * PhysConst::c / m_beta_boost;
287  return current_z_lab;
288  }
294  [[nodiscard]] amrex::Real dz_lab (amrex::Real dt, amrex::Real ref_ratio) const;
300  [[nodiscard]] int k_index_zlab (int i_buffer, int lev) const;
306  bool buffer_full (int i_buffer) {
307  return (k_index_zlab(i_buffer,0) == m_buffer_box[i_buffer].smallEnd(m_moving_window_dir));
308  }
309 
314  bool buffer_empty (int i_buffer) {
315  return ( m_buffer_counter[i_buffer] == 0) ;
316  }
317 
320 
324  void NullifyFirstFlush(int i_buffer) {m_first_flush_after_restart[i_buffer] = 0;}
325 
330 
334  void ResetBufferCounter(int i_buffer) {
335  m_buffer_counter[i_buffer] = 0;
336  }
340  void IncrementBufferFlushCounter(int i_buffer) {
341  m_buffer_flush_counter[i_buffer]++;
342  }
348  void SetSnapshotFullStatus (int i_buffer);
353 #ifdef WARPX_DIM_RZ
355  "Br", "Bt", "Bz",
356  "jr", "jt", "jz", "rho"};
358  "Br", "Bt", "Bz",
359  "jr", "jt", "jz",
360  "rho"};
361 #else
363  "Bx", "By", "Bz",
364  "jx", "jy", "jz", "rho"};
365 #endif
366 
367 
371  void MergeBuffersForPlotfile (int i_snapshot);
375  void InterleaveBufferAndSnapshotHeader (const std::string& buffer_Header,
376  const std::string& snapshot_Header);
380  void InterleaveFabArrayHeader (const std::string& Buffer_FabHeader_path,
381  const std::string& snapshot_FabHeader_path,
382  const std::string& newsnapshot_FabFilename);
386  void InterleaveSpeciesHeader(const std::string& buffer_species_Header_path,
387  const std::string& snapshot_species_Header_path,
388  const std::string& species_name, int new_data_index);
389 
393  void InterleaveParticleDataHeader(const std::string& buffer_ParticleHdrFilename,
394  const std::string& snapshot_ParticleHdrFilename);
397  void InitializeParticleFunctors () override;
398 
400  void ResetTotalParticlesInBuffer(int i_buffer);
402  void ClearParticleBuffer(int i_buffer);
404  void RedistributeParticleBuffer (int i_buffer);
406 
407  amrex::Real gettlab (int i_buffer) override {return m_t_lab[i_buffer];}
408  void settlab (int i_buffer, amrex::Real tlab) override {m_t_lab[i_buffer] = tlab; }
409  int get_buffer_k_index_hi (int i_buffer) override {return m_buffer_k_index_hi[i_buffer]; }
410  void set_buffer_k_index_hi (int i_buffer, int kindex) override {m_buffer_k_index_hi[i_buffer] = kindex;}
411  amrex::Real get_snapshot_domain_lo (int i_buffer, int idim) override {return m_snapshot_domain_lab[i_buffer].lo(idim); }
412  amrex::Real get_snapshot_domain_hi (int i_buffer, int idim) override {return m_snapshot_domain_lab[i_buffer].hi(idim); }
413  int get_flush_counter (int i_buffer) override {return m_buffer_flush_counter[i_buffer]; }
414  void set_flush_counter ([[maybe_unused]] int i_buffer, int flush_counter) override { m_buffer_flush_counter[i_buffer] = flush_counter; }
415  int get_last_valid_Zslice ( int i_buffer) override { return m_lastValidZSlice[i_buffer]; }
416  void set_last_valid_Zslice ( int i_buffer, int last_valid_Zslice) override { m_lastValidZSlice[i_buffer] = last_valid_Zslice; }
417  int get_snapshot_full_flag ( int i_buffer) override { return m_snapshot_full[i_buffer]; }
418  void set_snapshot_full ( int i_buffer, int snapshot_full) override { m_snapshot_full[i_buffer] = snapshot_full; }
419 };
420 #endif // WARPX_BTDIAGNOSTICS_H_
Definition: BTDiagnostics.H:28
amrex::Real m_beta_boost
Definition: BTDiagnostics.H:158
int m_moving_window_dir
Definition: BTDiagnostics.H:163
void DerivedInitData() override
Definition: BTDiagnostics.cpp:64
void InitializeParticleBuffer() override
Definition: BTDiagnostics.cpp:1432
amrex::Vector< std::string > m_cellcenter_varnames
Definition: BTDiagnostics.H:354
amrex::Vector< std::unique_ptr< amrex::MultiFab > > m_cell_centered_data
Definition: BTDiagnostics.H:244
void Flush(int i_buffer, bool force_flush) override
Flush m_mf_output and particles to file. Currently, a temporary customized output format for the buff...
Definition: BTDiagnostics.cpp:1006
void DefineCellCenteredMultiFab(int lev)
Definition: BTDiagnostics.cpp:505
void InitializeFieldFunctorsRZopenPMD(int lev) override
Definition: BTDiagnostics.cpp:656
amrex::Real dz_lab(amrex::Real dt, amrex::Real ref_ratio) const
Definition: BTDiagnostics.cpp:867
void IncrementBufferFlushCounter(int i_buffer)
Definition: BTDiagnostics.H:340
void InterleaveFabArrayHeader(const std::string &Buffer_FabHeader_path, const std::string &snapshot_FabHeader_path, const std::string &newsnapshot_FabFilename)
Definition: BTDiagnostics.cpp:1340
amrex::Real m_dz_snapshots_lab
Definition: BTDiagnostics.H:172
bool DoDump(int step, int i_buffer, bool force_flush=false) override
Definition: BTDiagnostics.cpp:294
bool GetKIndexInSnapshotBoxFlag(int i_buffer, int lev)
Definition: BTDiagnostics.cpp:999
amrex::Vector< std::string > m_cellcenter_varnames_fields
Definition: BTDiagnostics.H:357
void AddRZModesToOutputNames(const std::string &field, int ncomp, bool cellcenter_data)
Definition: BTDiagnostics.cpp:713
amrex::Real m_dt_snapshots_lab
Definition: BTDiagnostics.H:168
amrex::Vector< int > m_lastValidZSlice
Definition: BTDiagnostics.H:230
void set_last_valid_Zslice(int i_buffer, int last_valid_Zslice) override
Definition: BTDiagnostics.H:416
void PrepareFieldDataForOutput() override
Definition: BTDiagnostics.cpp:786
bool m_do_back_transformed_particles
Definition: BTDiagnostics.H:148
void DefineFieldBufferMultiFab(int i_buffer, int lev)
Definition: BTDiagnostics.cpp:898
amrex::Vector< amrex::RealBox > m_buffer_domain_lab
Definition: BTDiagnostics.H:181
void ReadParameters()
Definition: BTDiagnostics.cpp:206
amrex::Vector< int > m_buffer_counter
Definition: BTDiagnostics.H:216
amrex::Vector< amrex::Real > m_old_z_boost
Definition: BTDiagnostics.H:198
void set_snapshot_full(int i_buffer, int snapshot_full) override
Definition: BTDiagnostics.H:418
amrex::Vector< int > m_snapshot_full
Definition: BTDiagnostics.H:225
int get_snapshot_full_flag(int i_buffer) override
Definition: BTDiagnostics.H:417
void NullifyFirstFlush(int i_buffer)
Definition: BTDiagnostics.H:324
amrex::Vector< amrex::Vector< std::unique_ptr< ComputeDiagFunctor const > > > m_cell_center_functors
Definition: BTDiagnostics.H:248
amrex::Real gettlab(int i_buffer) override
Definition: BTDiagnostics.H:407
utils::parser::BTDIntervalsParser m_intervals
Definition: BTDiagnostics.H:41
void ClearParticleBuffer(int i_buffer)
Definition: BTDiagnostics.cpp:1500
void MergeBuffersForPlotfile(int i_snapshot)
Definition: BTDiagnostics.cpp:1124
void settlab(int i_buffer, amrex::Real tlab) override
Definition: BTDiagnostics.H:408
int m_buffer_size
Definition: BTDiagnostics.H:175
amrex::Vector< amrex::Real > m_current_z_lab
Definition: BTDiagnostics.H:192
amrex::Real get_snapshot_domain_hi(int i_buffer, int idim) override
Definition: BTDiagnostics.H:412
bool buffer_full(int i_buffer)
Definition: BTDiagnostics.H:306
amrex::Vector< int > m_field_buffer_multifab_defined
Definition: BTDiagnostics.H:329
void InterleaveBufferAndSnapshotHeader(const std::string &buffer_Header, const std::string &snapshot_Header)
Definition: BTDiagnostics.cpp:1295
amrex::Vector< int > m_max_buffer_multifabs
Definition: BTDiagnostics.H:221
bool buffer_empty(int i_buffer)
Definition: BTDiagnostics.H:314
void set_buffer_k_index_hi(int i_buffer, int kindex) override
Definition: BTDiagnostics.H:410
void InterleaveSpeciesHeader(const std::string &buffer_species_Header_path, const std::string &snapshot_species_Header_path, const std::string &species_name, int new_data_index)
Definition: BTDiagnostics.cpp:1370
void RedistributeParticleBuffer(int i_buffer)
Definition: BTDiagnostics.cpp:1117
void SetSnapshotFullStatus(int i_buffer)
Definition: BTDiagnostics.cpp:889
amrex::Vector< int > m_buffer_flush_counter
Definition: BTDiagnostics.H:233
bool DoComputeAndPack(int step, bool force_flush=false) override
Definition: BTDiagnostics.cpp:322
amrex::Vector< amrex::Box > m_buffer_box
Definition: BTDiagnostics.H:189
void InterleaveParticleDataHeader(const std::string &buffer_ParticleHdrFilename, const std::string &snapshot_ParticleHdrFilename)
Definition: BTDiagnostics.cpp:1394
bool m_plot_raw_fields_guards
Definition: BTDiagnostics.H:37
amrex::Real UpdateCurrentZBoostCoordinate(amrex::Real t_lab, amrex::Real t_boost) const
Definition: BTDiagnostics.H:275
void UpdateVarnamesForRZopenPMD()
Definition: BTDiagnostics.cpp:597
void PrepareParticleDataForOutput() override
Definition: BTDiagnostics.cpp:1453
amrex::Vector< amrex::Real > m_current_z_boost
Definition: BTDiagnostics.H:195
void UpdateBufferData() override
Definition: BTDiagnostics.cpp:764
void InitializeParticleFunctors() override
Definition: BTDiagnostics.cpp:1415
void InitializeBufferData(int i_buffer, int lev, bool restart=false) override
Definition: BTDiagnostics.cpp:333
void ResetTotalParticlesInBuffer(int i_buffer)
Definition: BTDiagnostics.cpp:1491
void DefineSnapshotGeometry(int i_buffer, int lev)
Definition: BTDiagnostics.cpp:958
BTDiagnostics(int i, const std::string &name)
Definition: BTDiagnostics.cpp:57
bool m_plot_raw_fields
Definition: BTDiagnostics.H:35
amrex::Vector< int > m_buffer_k_index_hi
Definition: BTDiagnostics.H:237
bool m_do_back_transformed_fields
Definition: BTDiagnostics.H:144
amrex::Real m_gamma_boost
Definition: BTDiagnostics.H:157
amrex::Vector< amrex::IntVect > m_snapshot_ncells_lab
Definition: BTDiagnostics.H:183
amrex::Real UpdateCurrentZLabCoordinate(amrex::Real t_lab, amrex::Real t_boost) const
Definition: BTDiagnostics.H:284
int k_index_zlab(int i_buffer, int lev) const
Definition: BTDiagnostics.cpp:874
amrex::Vector< amrex::Vector< amrex::Geometry > > m_geom_snapshot
Definition: BTDiagnostics.H:209
int get_last_valid_Zslice(int i_buffer) override
Definition: BTDiagnostics.H:415
amrex::Real get_snapshot_domain_lo(int i_buffer, int idim) override
Definition: BTDiagnostics.H:411
int get_flush_counter(int i_buffer) override
Definition: BTDiagnostics.H:413
void ResetBufferCounter(int i_buffer)
Definition: BTDiagnostics.H:334
amrex::Vector< int > m_first_flush_after_restart
Definition: BTDiagnostics.H:319
amrex::Vector< amrex::Box > m_snapshot_box
Definition: BTDiagnostics.H:186
void PrepareBufferData() override
Definition: BTDiagnostics.cpp:741
void InitializeFieldFunctors(int lev) override
Definition: BTDiagnostics.cpp:526
int m_num_snapshots_lab
Definition: BTDiagnostics.H:166
amrex::Vector< amrex::Real > m_t_lab
Definition: BTDiagnostics.H:178
void set_flush_counter([[maybe_unused]] int i_buffer, int flush_counter) override
Definition: BTDiagnostics.H:414
int get_buffer_k_index_hi(int i_buffer) override
Definition: BTDiagnostics.H:409
amrex::Vector< int > m_snapshot_geometry_defined
Definition: BTDiagnostics.H:327
bool GetZSliceInDomainFlag(int i_buffer, int lev)
Definition: BTDiagnostics.cpp:981
base class for diagnostics. Contains main routines to filter, compute and flush diagnostics.
Definition: Diagnostics.H:31
amrex::Vector< amrex::RealBox > m_snapshot_domain_lab
Definition: Diagnostics.H:319
This class is a parser for multiple slices of the form x,y,z,... where x, y and z are slices of the f...
Definition: IntervalsParser.H:177
static constexpr auto c
vacuum speed of light [m/s]
Definition: constant.H:44
i
Definition: check_interp_points_and_weights.py:174
name
Definition: run_automated.py:229
float dt
Definition: stencil.py:442
string field
Definition: video_yt.py:31