aboutsummaryrefslogtreecommitdiff
path: root/Source/Diagnostics/BoostedFrameDiagnostic.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'Source/Diagnostics/BoostedFrameDiagnostic.cpp')
-rw-r--r--Source/Diagnostics/BoostedFrameDiagnostic.cpp729
1 files changed, 551 insertions, 178 deletions
diff --git a/Source/Diagnostics/BoostedFrameDiagnostic.cpp b/Source/Diagnostics/BoostedFrameDiagnostic.cpp
index 45343a0cb..297b4f5be 100644
--- a/Source/Diagnostics/BoostedFrameDiagnostic.cpp
+++ b/Source/Diagnostics/BoostedFrameDiagnostic.cpp
@@ -1,6 +1,8 @@
#include <AMReX_MultiFabUtil.H>
+#include <AMReX_MultiFabUtil_C.H>
#include "BoostedFrameDiagnostic.H"
+#include "SliceDiagnostic.H"
#include "WarpX_f.H"
#include "WarpX.H"
@@ -16,7 +18,9 @@ using namespace amrex;
*/
namespace
{
- const std::vector<std::string> particle_field_names = {"w", "x", "y", "z", "ux", "uy", "uz"};
+
+ const std::vector<std::string> particle_field_names = {"w", "x", "y",
+ "z", "ux", "uy", "uz"};
/*
Creates the HDF5 file in truncate mode and closes it.
@@ -336,7 +340,8 @@ namespace
*/
void output_write_field(const std::string& file_path,
const std::string& field_path,
- const MultiFab& mf, const int comp)
+ const MultiFab& mf, const int comp,
+ const int lo_x, const int lo_y, const int lo_z)
{
BL_PROFILE("output_write_field");
@@ -374,8 +379,15 @@ namespace
// slab lo index and shape.
#if (AMREX_SPACEDIM == 3)
hsize_t slab_offsets[3], slab_dims[3];
+ int shift[3];
+ shift[0] = lo_x;
+ shift[1] = lo_y;
+ shift[2] = lo_z;
#else
hsize_t slab_offsets[2], slab_dims[2];
+ int shift[2];
+ shift[0] = lo_x;
+ shift[1] = lo_z;
#endif
hid_t slab_dataspace;
@@ -396,7 +408,7 @@ namespace
{
AMREX_ASSERT(lo_vec[idim] >= 0);
AMREX_ASSERT(hi_vec[idim] > lo_vec[idim]);
- slab_offsets[idim] = lo_vec[idim];
+ slab_offsets[idim] = lo_vec[idim] - shift[idim];
slab_dims[idim] = hi_vec[idim] - lo_vec[idim] + 1;
}
@@ -444,39 +456,14 @@ namespace
}
#endif
-namespace
+bool compare_tlab_uptr(const std::unique_ptr<LabFrameDiag>&a,
+ const std::unique_ptr<LabFrameDiag>&b)
{
- void
- CopySlice(MultiFab& tmp, MultiFab& buf, int k_lab,
- const Gpu::ManagedDeviceVector<int>& map_actual_fields_to_dump)
- {
- const int ncomp_to_dump = map_actual_fields_to_dump.size();
- // Copy data from MultiFab tmp to MultiFab data_buffer[i].
-#ifdef _OPENMP
-#pragma omp parallel if (Gpu::notInLaunchRegion())
-#endif
- for (MFIter mfi(tmp, TilingIfNotGPU()); mfi.isValid(); ++mfi) {
- Array4< Real> tmp_arr = tmp[mfi].array();
- Array4< Real> buf_arr = buf[mfi].array();
- // For 3D runs, tmp is a 2D (x,y) multifab, that contains only
- // slice to write to file.
- const Box& bx = mfi.tilebox();
-
- const auto field_map_ptr = map_actual_fields_to_dump.dataPtr();
- ParallelFor(bx, ncomp_to_dump,
- [=] AMREX_GPU_DEVICE (int i, int j, int k, int n)
- {
- const int icomp = field_map_ptr[n];
-#if (AMREX_SPACEDIM == 3)
- buf_arr(i,j,k_lab,n) = tmp_arr(i,j,k,icomp);
-#else
- buf_arr(i,k_lab,k,n) = tmp_arr(i,j,k,icomp);
-#endif
- }
- );
- }
- }
+ return a->t_lab < b->t_lab;
+}
+namespace
+{
void
LorentzTransformZ(MultiFab& data, Real gamma_boost, Real beta_boost, int ncomp)
{
@@ -497,21 +484,27 @@ LorentzTransformZ(MultiFab& data, Real gamma_boost, Real beta_boost, int ncomp)
// Transform the transverse E and B fields. Note that ez and bz are not
// changed by the tranform.
Real e_lab, b_lab, j_lab, r_lab;
- e_lab = gamma_boost * (arr(i, j, k, 0) + beta_boost*clight*arr(i, j, k, 4));
- b_lab = gamma_boost * (arr(i, j, k, 4) + beta_boost*arr(i, j, k, 0)/clight);
+ e_lab = gamma_boost * (arr(i, j, k, 0) +
+ beta_boost*clight*arr(i, j, k, 4));
+ b_lab = gamma_boost * (arr(i, j, k, 4) +
+ beta_boost*arr(i, j, k, 0)/clight);
arr(i, j, k, 0) = e_lab;
arr(i, j, k, 4) = b_lab;
- e_lab = gamma_boost * (arr(i, j, k, 1) - beta_boost*clight*arr(i, j, k, 3));
- b_lab = gamma_boost * (arr(i, j, k, 3) - beta_boost*arr(i, j, k, 1)/clight);
+ e_lab = gamma_boost * (arr(i, j, k, 1) -
+ beta_boost*clight*arr(i, j, k, 3));
+ b_lab = gamma_boost * (arr(i, j, k, 3) -
+ beta_boost*arr(i, j, k, 1)/clight);
arr(i, j, k, 1) = e_lab;
arr(i, j, k, 3) = b_lab;
// Transform the charge and current density. Only the z component of j is affected.
- j_lab = gamma_boost*(arr(i, j, k, 8) + beta_boost*clight*arr(i, j, k, 9));
- r_lab = gamma_boost*(arr(i, j, k, 9) + beta_boost*arr(i, j, k, 8)/clight);
+ j_lab = gamma_boost*(arr(i, j, k, 8) +
+ beta_boost*clight*arr(i, j, k, 9));
+ r_lab = gamma_boost*(arr(i, j, k, 9) +
+ beta_boost*arr(i, j, k, 8)/clight);
arr(i, j, k, 8) = j_lab;
arr(i, j, k, 9) = r_lab;
@@ -524,15 +517,20 @@ LorentzTransformZ(MultiFab& data, Real gamma_boost, Real beta_boost, int ncomp)
BoostedFrameDiagnostic::
BoostedFrameDiagnostic(Real zmin_lab, Real zmax_lab, Real v_window_lab,
Real dt_snapshots_lab, int N_snapshots,
+ Real dt_slice_snapshots_lab, int N_slice_snapshots,
Real gamma_boost, Real t_boost, Real dt_boost,
- int boost_direction, const Geometry& geom)
+ int boost_direction, const Geometry& geom,
+ amrex::RealBox& slice_realbox)
: gamma_boost_(gamma_boost),
dt_snapshots_lab_(dt_snapshots_lab),
dt_boost_(dt_boost),
N_snapshots_(N_snapshots),
- boost_direction_(boost_direction)
+ boost_direction_(boost_direction),
+ N_slice_snapshots_(N_slice_snapshots),
+ dt_slice_snapshots_lab_(dt_slice_snapshots_lab)
{
+
BL_PROFILE("BoostedFrameDiagnostic::BoostedFrameDiagnostic");
AMREX_ALWAYS_ASSERT(WarpX::do_boosted_frame_fields or
@@ -553,12 +551,8 @@ BoostedFrameDiagnostic(Real zmin_lab, Real zmax_lab, Real v_window_lab,
// Ny_lab = 1;
IntVect prob_ncells_lab = {Nx_lab, Nz_lab};
#endif
-
writeMetaData();
- if (WarpX::do_boosted_frame_fields) data_buffer_.resize(N_snapshots);
- if (WarpX::do_boosted_frame_particles) particles_buffer_.resize(N_snapshots);
-
// Query fields to dump
std::vector<std::string> user_fields_to_dump;
ParmParse pp("warpx");
@@ -581,6 +575,9 @@ BoostedFrameDiagnostic(Real zmin_lab, Real zmax_lab, Real v_window_lab,
}
}
+ // allocating array with total number of lab frame diags (snapshots+slices)
+ LabFrameDiags_.resize(N_snapshots+N_slice_snapshots);
+
for (int i = 0; i < N_snapshots; ++i) {
Real t_lab = i * dt_snapshots_lab_;
// Get simulation domain physical coordinates (in boosted frame).
@@ -589,14 +586,95 @@ BoostedFrameDiagnostic(Real zmin_lab, Real zmax_lab, Real v_window_lab,
// x and y bounds are the same for lab frame and boosted frame
prob_domain_lab.setLo(AMREX_SPACEDIM-1, zmin_lab + v_window_lab * t_lab);
prob_domain_lab.setHi(AMREX_SPACEDIM-1, zmax_lab + v_window_lab * t_lab);
- // Construct LabSnapShot
- LabSnapShot snapshot(t_lab, t_boost, prob_domain_lab,
- prob_ncells_lab, ncomp_to_dump,
- mesh_field_names, i, *this);
- snapshots_.push_back(snapshot);
- buff_counter_.push_back(0);
- if (WarpX::do_boosted_frame_fields) data_buffer_[i].reset( nullptr );
+ Box diag_box = geom.Domain();
+ LabFrameDiags_[i].reset(new LabFrameSnapShot(t_lab, t_boost,
+ inv_gamma_boost_, inv_beta_boost_, dz_lab_,
+ prob_domain_lab, prob_ncells_lab,
+ ncomp_to_dump, mesh_field_names, prob_domain_lab,
+ diag_box, i));
+ }
+
+
+ for (int i = 0; i < N_slice_snapshots; ++i) {
+
+ amrex::Real cell_dx = 0;
+ amrex::Real cell_dy = 0;
+ IntVect slice_ncells_lab ;
+
+ // To construct LabFrameSlice(), the location of lo() and hi() of the
+ // reduced diag is computed using the user-defined values of the
+ // reduced diag (1D, 2D, or 3D). For visualization of the diagnostics,
+ // the number of cells in each dimension is required and
+ // is computed below for the reduced back-transformed lab-frame diag,
+ // similar to the full-diag.
+ const amrex::Real* current_slice_lo = slice_realbox.lo();
+ const amrex::Real* current_slice_hi = slice_realbox.hi();
+
+ const amrex::Real zmin_slice_lab = current_slice_lo[AMREX_SPACEDIM-1] /
+ ( (1.+beta_boost_)*gamma_boost_);
+ const amrex::Real zmax_slice_lab = current_slice_hi[AMREX_SPACEDIM-1] /
+ ( (1.+beta_boost_)*gamma_boost_);
+ int Nz_slice_lab = static_cast<unsigned>((
+ zmax_slice_lab - zmin_slice_lab) * inv_dz_lab_);
+ int Nx_slice_lab = ( current_slice_hi[0] - current_slice_lo[0] ) /
+ geom.CellSize(0);
+ if (Nx_slice_lab == 0 ) Nx_slice_lab = 1;
+ // if the x-dimension is reduced, increase total_cells by 1
+ // to be consistent with the number of cells created for the output.
+ if (Nx_lab != Nx_slice_lab) Nx_slice_lab++;
+ cell_dx = geom.CellSize(0);
+#if (AMREX_SPACEDIM == 3)
+ int Ny_slice_lab = ( current_slice_hi[1] - current_slice_lo[1]) /
+ geom.CellSize(1);
+ if (Ny_slice_lab == 0 ) Ny_slice_lab = 1;
+ // if the y-dimension is reduced, increase total_cells by 1
+ // to be consistent with the number of cells created for the output.
+ if (Ny_lab != Ny_slice_lab) Ny_slice_lab++;
+ slice_ncells_lab = {Nx_slice_lab, Ny_slice_lab, Nz_slice_lab};
+ cell_dy = geom.CellSize(1);
+#else
+ slice_ncells_lab = {Nx_slice_lab, Nz_slice_lab};
+#endif
+
+ IntVect slice_lo(AMREX_D_DECL(0,0,0));
+ IntVect slice_hi(AMREX_D_DECL(1,1,1));
+
+ for ( int i_dim=0; i_dim<AMREX_SPACEDIM; ++i_dim)
+ {
+ slice_lo[i_dim] = (slice_realbox.lo(i_dim) - geom.ProbLo(i_dim) -
+ 0.5*geom.CellSize(i_dim))/geom.CellSize(i_dim);
+ slice_hi[i_dim] = (slice_realbox.hi(i_dim) - geom.ProbLo(i_dim) -
+ 0.5*geom.CellSize(i_dim))/geom.CellSize(i_dim);
+ if (slice_lo[i_dim] == slice_hi[i_dim])
+ {
+ slice_hi[i_dim] = slice_lo[i_dim] + 1;
+ }
+ }
+ Box stmp(slice_lo,slice_hi);
+ Box slicediag_box = stmp;
+
+ Real t_slice_lab = i * dt_slice_snapshots_lab_ ;
+ RealBox prob_domain_lab = geom.ProbDomain();
+ // replace z bounds by lab-frame coordinates
+ prob_domain_lab.setLo(AMREX_SPACEDIM-1, zmin_lab + v_window_lab * t_slice_lab);
+ prob_domain_lab.setHi(AMREX_SPACEDIM-1, zmax_lab + v_window_lab * t_slice_lab);
+ RealBox slice_dom_lab = slice_realbox;
+ // replace z bounds of slice in lab-frame coordinates
+ // note : x and y bounds are the same for lab and boosted frames
+ // initial lab slice extent //
+ slice_dom_lab.setLo(AMREX_SPACEDIM-1, zmin_slice_lab + v_window_lab * t_slice_lab );
+ slice_dom_lab.setHi(AMREX_SPACEDIM-1, zmax_slice_lab +
+ v_window_lab * t_slice_lab );
+
+ // construct labframeslice
+ LabFrameDiags_[i+N_snapshots].reset(new LabFrameSlice(t_slice_lab, t_boost,
+ inv_gamma_boost_, inv_beta_boost_, dz_lab_,
+ prob_domain_lab, slice_ncells_lab,
+ ncomp_to_dump, mesh_field_names, slice_dom_lab,
+ slicediag_box, i, cell_dx, cell_dy));
}
+ // sort diags based on their respective t_lab
+ std::stable_sort(LabFrameDiags_.begin(), LabFrameDiags_.end(), compare_tlab_uptr);
AMREX_ALWAYS_ASSERT(max_box_size_ >= num_buffer_);
}
@@ -612,18 +690,19 @@ void BoostedFrameDiagnostic::Flush(const Geometry& geom)
const std::vector<std::string> species_names = mypc.GetSpeciesNames();
// Loop over BFD snapshots
- for (int i = 0; i < N_snapshots_; ++i) {
+ for (int i = 0; i < LabFrameDiags_.size(); ++i) {
- Real zmin_lab = snapshots_[i].prob_domain_lab_.lo(AMREX_SPACEDIM-1);
- int i_lab = (snapshots_[i].current_z_lab - zmin_lab) / dz_lab_;
+ Real zmin_lab = LabFrameDiags_[i]->prob_domain_lab_.lo(AMREX_SPACEDIM-1);
+ int i_lab = (LabFrameDiags_[i]->current_z_lab - zmin_lab) / dz_lab_;
- if (buff_counter_[i] != 0) {
+ if (LabFrameDiags_[i]->buff_counter_ != 0) {
if (WarpX::do_boosted_frame_fields) {
- const BoxArray& ba = data_buffer_[i]->boxArray();
+ const BoxArray& ba = LabFrameDiags_[i]->data_buffer_->boxArray();
const int hi = ba[0].bigEnd(boost_direction_);
- const int lo = hi - buff_counter_[i] + 1;
+ const int lo = hi - LabFrameDiags_[i]->buff_counter_ + 1;
- Box buff_box = geom.Domain();
+ //Box buff_box = geom.Domain();
+ Box buff_box = LabFrameDiags_[i]->buff_box_;
buff_box.setSmall(boost_direction_, lo);
buff_box.setBig(boost_direction_, hi);
@@ -631,18 +710,23 @@ void BoostedFrameDiagnostic::Flush(const Geometry& geom)
buff_ba.maxSize(max_box_size_);
DistributionMapping buff_dm(buff_ba);
- const int ncomp = data_buffer_[i]->nComp();
+ const int ncomp = LabFrameDiags_[i]->data_buffer_->nComp();
MultiFab tmp(buff_ba, buff_dm, ncomp, 0);
- tmp.copy(*data_buffer_[i], 0, 0, ncomp);
+ tmp.copy(*LabFrameDiags_[i]->data_buffer_, 0, 0, ncomp);
#ifdef WARPX_USE_HDF5
- for (int comp = 0; comp < ncomp; ++comp)
- output_write_field(snapshots_[i].file_name, mesh_field_names[comp], tmp, comp);
+ for (int comp = 0; comp < ncomp; ++comp) {
+ output_write_field(LabFrameDiags_[i]->file_name,
+ mesh_field_names[comp], tmp, comp,
+ lbound(buff_box).x, lbound(buff_box).y,
+ lbound(buff_box).z);
+ }
#else
std::stringstream ss;
- ss << snapshots_[i].file_name << "/Level_0/" << Concatenate("buffer", i_lab, 5);
+ ss << LabFrameDiags_[i]->file_name << "/Level_0/"
+ << Concatenate("buffer", i_lab, 5);
VisMF::Write(tmp, ss.str());
#endif
}
@@ -655,25 +739,31 @@ void BoostedFrameDiagnostic::Flush(const Geometry& geom)
species_names[mypc.mapSpeciesBoostedFrameDiags(j)];
#ifdef WARPX_USE_HDF5
// Dump species data
- writeParticleDataHDF5(particles_buffer_[i][j],
- snapshots_[i].file_name,
+ writeParticleDataHDF5(LabFrameDiags_[i]->particles_buffer_[j],
+ LabFrameDiags_[i]->file_name,
species_name);
#else
std::stringstream part_ss;
- part_ss << snapshots_[i].file_name + "/" + species_name + "/";
+ part_ss << LabFrameDiags_[i]->file_name + "/" +
+ species_name + "/";
// Dump species data
- writeParticleData(particles_buffer_[i][j], part_ss.str(), i_lab);
+ writeParticleData(LabFrameDiags_[i]->particles_buffer_[j],
+ part_ss.str(), i_lab);
#endif
}
- particles_buffer_[i].clear();
+ LabFrameDiags_[i]->particles_buffer_.clear();
}
- buff_counter_[i] = 0;
+ LabFrameDiags_[i]->buff_counter_ = 0;
}
}
VisMF::SetHeaderVersion(current_version);
}
+
+
+
+
void
BoostedFrameDiagnostic::
writeLabFrameData(const MultiFab* cell_centered_data,
@@ -681,7 +771,6 @@ writeLabFrameData(const MultiFab* cell_centered_data,
const Geometry& geom, const Real t_boost, const Real dt) {
BL_PROFILE("BoostedFrameDiagnostic::writeLabFrameData");
-
VisMF::Header::Version current_version = VisMF::GetHeaderVersion();
VisMF::SetHeaderVersion(amrex::VisMF::Header::NoFabHeader_v1);
@@ -690,97 +779,140 @@ writeLabFrameData(const MultiFab* cell_centered_data,
const Real zhi_boost = domain_z_boost.hi(boost_direction_);
const std::vector<std::string> species_names = mypc.GetSpeciesNames();
+ Real prev_t_lab = -dt;
+ std::unique_ptr<amrex::MultiFab> tmp_slice_ptr;
+ std::unique_ptr<amrex::MultiFab> slice;
+ amrex::Vector<WarpXParticleContainer::DiagnosticParticleData> tmp_particle_buffer;
// Loop over snapshots
- for (int i = 0; i < N_snapshots_; ++i) {
-
+ for (int i = 0; i < LabFrameDiags_.size(); ++i) {
// Get updated z position of snapshot
- const Real old_z_boost = snapshots_[i].current_z_boost;
- snapshots_[i].updateCurrentZPositions(t_boost,
+ const Real old_z_boost = LabFrameDiags_[i]->current_z_boost;
+ LabFrameDiags_[i]->updateCurrentZPositions(t_boost,
inv_gamma_boost_,
inv_beta_boost_);
- Real zmin_lab = snapshots_[i].prob_domain_lab_.lo(AMREX_SPACEDIM-1);
- Real zmax_lab = snapshots_[i].prob_domain_lab_.hi(AMREX_SPACEDIM-1);
+ Real diag_zmin_lab = LabFrameDiags_[i]->diag_domain_lab_.lo(AMREX_SPACEDIM-1);
+ Real diag_zmax_lab = LabFrameDiags_[i]->diag_domain_lab_.hi(AMREX_SPACEDIM-1);
- // If snapshot out of the domain, nothing to do
- if ( (snapshots_[i].current_z_boost < zlo_boost) or
- (snapshots_[i].current_z_boost > zhi_boost) or
- (snapshots_[i].current_z_lab < zmin_lab) or
- (snapshots_[i].current_z_lab > zmax_lab) ) continue;
+ if ( ( LabFrameDiags_[i]->current_z_boost < zlo_boost) or
+ ( LabFrameDiags_[i]->current_z_boost > zhi_boost) or
+ ( LabFrameDiags_[i]->current_z_lab < diag_zmin_lab) or
+ ( LabFrameDiags_[i]->current_z_lab > diag_zmax_lab) ) continue;
// Get z index of data_buffer_ (i.e. in the lab frame) where
// simulation domain (t', [zmin',zmax']), back-transformed to lab
// frame, intersects with snapshot.
- int i_lab = (snapshots_[i].current_z_lab - zmin_lab) / dz_lab_;
-
+ Real dom_zmin_lab = LabFrameDiags_[i]->prob_domain_lab_.lo(AMREX_SPACEDIM-1);
+ int i_lab = ( LabFrameDiags_[i]->current_z_lab - dom_zmin_lab) / dz_lab_;
// If buffer of snapshot i is empty...
- if (buff_counter_[i] == 0) {
- // ... reset fields buffer data_buffer_[i]
+ if ( LabFrameDiags_[i]->buff_counter_ == 0) {
+ // ... reset fields buffer data_buffer_
if (WarpX::do_boosted_frame_fields) {
- Box buff_box = geom.Domain();
- buff_box.setSmall(boost_direction_, i_lab - num_buffer_ + 1);
- buff_box.setBig(boost_direction_, i_lab);
- BoxArray buff_ba(buff_box);
+ LabFrameDiags_[i]->buff_box_.setSmall(boost_direction_,
+ i_lab - num_buffer_ + 1);
+ LabFrameDiags_[i]->buff_box_.setBig(boost_direction_, i_lab);
+
+ BoxArray buff_ba(LabFrameDiags_[i]->buff_box_);
buff_ba.maxSize(max_box_size_);
DistributionMapping buff_dm(buff_ba);
- data_buffer_[i].reset( new MultiFab(buff_ba, buff_dm, ncomp_to_dump, 0) );
+ LabFrameDiags_[i]->data_buffer_.reset( new MultiFab(buff_ba,
+ buff_dm, ncomp_to_dump, 0) );
}
// ... reset particle buffer particles_buffer_[i]
if (WarpX::do_boosted_frame_particles)
- particles_buffer_[i].resize(mypc.nSpeciesBoostedFrameDiags());
+ LabFrameDiags_[i]->particles_buffer_.resize(
+ mypc.nSpeciesBoostedFrameDiags());
}
if (WarpX::do_boosted_frame_fields) {
const int ncomp = cell_centered_data->nComp();
const int start_comp = 0;
const bool interpolate = true;
- // Get slice in the boosted frame
- std::unique_ptr<MultiFab> slice = amrex::get_slice_data(boost_direction_,
- snapshots_[i].current_z_boost,
- *cell_centered_data, geom,
- start_comp, ncomp, interpolate);
-
- // transform it to the lab frame
- LorentzTransformZ(*slice, gamma_boost_, beta_boost_, ncomp);
- // Create a 2D box for the slice in the boosted frame
- Real dx = geom.CellSize(boost_direction_);
- int i_boost = (snapshots_[i].current_z_boost - geom.ProbLo(boost_direction_))/dx;
- Box slice_box = geom.Domain();
- slice_box.setSmall(boost_direction_, i_boost);
- slice_box.setBig(boost_direction_, i_boost);
- // Make it a BoxArray slice_ba
- BoxArray slice_ba(slice_box);
- slice_ba.maxSize(max_box_size_);
- // Create MultiFab tmp on slice_ba with data from slice
- MultiFab tmp(slice_ba, data_buffer_[i]->DistributionMap(), ncomp, 0);
- tmp.copy(*slice, 0, 0, ncomp);
-
- // Copy data from MultiFab tmp to MultiDab data_buffer[i]
- CopySlice(tmp, *data_buffer_[i], i_lab, map_actual_fields_to_dump);
+ // slice containing back-transformed data is generated only if t_lab != prev_t_lab and is re-used if multiple diags have the same z_lab,t_lab.
+ if (LabFrameDiags_[i]->t_lab != prev_t_lab ) {
+ if (slice)
+ {
+ slice.reset(new MultiFab);
+ slice.reset(nullptr);
+ }
+ slice = amrex::get_slice_data(boost_direction_,
+ LabFrameDiags_[i]->current_z_boost,
+ *cell_centered_data, geom,
+ start_comp, ncomp,
+ interpolate);
+ // Back-transform data to the lab-frame
+ LorentzTransformZ(*slice, gamma_boost_, beta_boost_, ncomp);
+ }
+ // Create a 2D box for the slice in the boosted frame
+ Real dx = geom.CellSize(boost_direction_);
+ int i_boost = ( LabFrameDiags_[i]->current_z_boost -
+ geom.ProbLo(boost_direction_))/dx;
+ //Box slice_box = geom.Domain();
+ Box slice_box = LabFrameDiags_[i]->buff_box_;
+ slice_box.setSmall(boost_direction_, i_boost);
+ slice_box.setBig(boost_direction_, i_boost);
+
+ // Make it a BoxArray slice_ba
+ BoxArray slice_ba(slice_box);
+ slice_ba.maxSize(max_box_size_);
+ tmp_slice_ptr = std::unique_ptr<MultiFab>(new MultiFab(slice_ba,
+ LabFrameDiags_[i]->data_buffer_->DistributionMap(),
+ ncomp, 0));
+
+ // slice is re-used if the t_lab of a diag is equal to
+ // that of the previous diag.
+ // Back-transformed data is copied from slice
+ // which has the dmap of the domain to
+ // tmp_slice_ptr which has the dmap of the
+ // data_buffer that stores the back-transformed data.
+ tmp_slice_ptr->copy(*slice, 0, 0, ncomp);
+ LabFrameDiags_[i]->AddDataToBuffer(*tmp_slice_ptr, i_lab,
+ map_actual_fields_to_dump);
+ tmp_slice_ptr.reset(new MultiFab);
+ tmp_slice_ptr.reset(nullptr);
}
if (WarpX::do_boosted_frame_particles) {
- mypc.GetLabFrameData(snapshots_[i].file_name, i_lab, boost_direction_,
- old_z_boost, snapshots_[i].current_z_boost,
- t_boost, snapshots_[i].t_lab, dt, particles_buffer_[i]);
- }
+ if (LabFrameDiags_[i]->t_lab != prev_t_lab ) {
+ if (tmp_particle_buffer.size()>0)
+ {
+ tmp_particle_buffer.clear();
+ tmp_particle_buffer.shrink_to_fit();
+ }
+ tmp_particle_buffer.resize(mypc.nSpeciesBoostedFrameDiags());
+ mypc.GetLabFrameData( LabFrameDiags_[i]->file_name, i_lab,
+ boost_direction_, old_z_boost,
+ LabFrameDiags_[i]->current_z_boost,
+ t_boost, LabFrameDiags_[i]->t_lab, dt,
+ tmp_particle_buffer);
+ }
+ LabFrameDiags_[i]->AddPartDataToParticleBuffer(tmp_particle_buffer,
+ mypc.nSpeciesBoostedFrameDiags());
+ }
- ++buff_counter_[i];
+ ++LabFrameDiags_[i]->buff_counter_;
+ prev_t_lab = LabFrameDiags_[i]->t_lab;
// If buffer full, write to disk.
- if (buff_counter_[i] == num_buffer_) {
+ if ( LabFrameDiags_[i]->buff_counter_ == num_buffer_) {
if (WarpX::do_boosted_frame_fields) {
#ifdef WARPX_USE_HDF5
- for (int comp = 0; comp < data_buffer_[i]->nComp(); ++comp)
- output_write_field(snapshots_[i].file_name, mesh_field_names[comp],
- *data_buffer_[i], comp);
+
+ Box buff_box = LabFrameDiags_[i]->buff_box_;
+ for (int comp = 0; comp < LabFrameDiags_[i]->data_buffer_->nComp(); ++comp)
+ output_write_field( LabFrameDiags_[i]->file_name,
+ mesh_field_names[comp],
+ *LabFrameDiags_[i]->data_buffer_, comp,
+ lbound(buff_box).x, lbound(buff_box).y,
+ lbound(buff_box).z);
#else
std::stringstream mesh_ss;
- mesh_ss << snapshots_[i].file_name << "/Level_0/" << Concatenate("buffer", i_lab, 5);
- VisMF::Write(*data_buffer_[i], mesh_ss.str());
+ mesh_ss << LabFrameDiags_[i]->file_name << "/Level_0/" <<
+ Concatenate("buffer", i_lab, 5);
+ VisMF::Write( (*LabFrameDiags_[i]->data_buffer_), mesh_ss.str());
#endif
}
@@ -788,24 +920,27 @@ writeLabFrameData(const MultiFab* cell_centered_data,
// Loop over species to be dumped to BFD
for (int j = 0; j < mypc.nSpeciesBoostedFrameDiags(); ++j) {
// Get species name
- const std::string species_name = species_names[mypc.mapSpeciesBoostedFrameDiags(j)];
+ const std::string species_name = species_names[
+ mypc.mapSpeciesBoostedFrameDiags(j)];
#ifdef WARPX_USE_HDF5
// Write data to disk (HDF5)
- writeParticleDataHDF5(particles_buffer_[i][j],
- snapshots_[i].file_name,
+ writeParticleDataHDF5(LabFrameDiags_[i]->particles_buffer_[j],
+ LabFrameDiags_[i]->file_name,
species_name);
#else
std::stringstream part_ss;
- part_ss << snapshots_[i].file_name + "/" + species_name + "/";
+ part_ss << LabFrameDiags_[i]->file_name + "/" +
+ species_name + "/";
// Write data to disk (custom)
- writeParticleData(particles_buffer_[i][j], part_ss.str(), i_lab);
+ writeParticleData(LabFrameDiags_[i]->particles_buffer_[j],
+ part_ss.str(), i_lab);
#endif
}
- particles_buffer_[i].clear();
+ LabFrameDiags_[i]->particles_buffer_.clear();
}
- buff_counter_[i] = 0;
+ LabFrameDiags_[i]->buff_counter_ = 0;
}
}
@@ -823,7 +958,8 @@ writeParticleDataHDF5(const WarpXParticleContainer::DiagnosticParticleData& pdat
Vector<long> particle_counts(ParallelDescriptor::NProcs(), 0);
Vector<long> particle_offsets(ParallelDescriptor::NProcs(), 0);
- ParallelAllGather::AllGather(np, particle_counts.data(), ParallelContext::CommunicatorAll());
+ ParallelAllGather::AllGather(np, particle_counts.data(),
+ ParallelContext::CommunicatorAll());
long total_np = 0;
for (int i = 0; i < ParallelDescriptor::NProcs(); ++i) {
@@ -854,7 +990,8 @@ writeParticleDataHDF5(const WarpXParticleContainer::DiagnosticParticleData& pdat
output_write_particle_field(name, field_path,
pdata.GetRealData(k).data(),
particle_counts[ParallelDescriptor::MyProc()],
- particle_offsets[ParallelDescriptor::MyProc()] + old_np);
+ particle_offsets[ParallelDescriptor::MyProc()]
+ + old_np);
}
}
#endif
@@ -871,7 +1008,6 @@ writeParticleData(const WarpXParticleContainer::DiagnosticParticleData& pdata,
const int MyProc = ParallelDescriptor::MyProc();
auto np = pdata.GetRealData(DiagIdx::w).size();
-
if (np == 0) return;
field_name = name + Concatenate("w_", i_lab, 5) + "_" + std::to_string(MyProc);
@@ -917,14 +1053,14 @@ writeMetaData ()
BL_PROFILE("BoostedFrameDiagnostic::writeMetaData");
if (ParallelDescriptor::IOProcessor()) {
-
- if (!UtilCreateDirectory(WarpX::lab_data_directory, 0755))
- CreateDirectoryFailed(WarpX::lab_data_directory);
+ const std::string fullpath = WarpX::lab_data_directory + "/snapshots";
+ if (!UtilCreateDirectory(fullpath, 0755))
+ CreateDirectoryFailed(fullpath);
VisMF::IO_Buffer io_buffer(VisMF::IO_Buffer_Size);
std::ofstream HeaderFile;
HeaderFile.rdbuf()->pubsetbuf(io_buffer.dataPtr(), io_buffer.size());
- std::string HeaderFileName(WarpX::lab_data_directory + "/Header");
+ std::string HeaderFileName(WarpX::lab_data_directory + "/snapshots/Header");
HeaderFile.open(HeaderFileName.c_str(), std::ofstream::out |
std::ofstream::trunc |
std::ofstream::binary);
@@ -937,31 +1073,81 @@ writeMetaData ()
HeaderFile << dt_snapshots_lab_ << "\n";
HeaderFile << gamma_boost_ << "\n";
HeaderFile << beta_boost_ << "\n";
+
+ if (N_slice_snapshots_ > 0) {
+ const std::string fullpath_slice = WarpX::lab_data_directory + "/slices";
+ if (!UtilCreateDirectory(fullpath_slice, 0755))
+ CreateDirectoryFailed(fullpath_slice);
+
+ VisMF::IO_Buffer io_buffer_slice(VisMF::IO_Buffer_Size);
+ std::ofstream HeaderFile_slice;
+ HeaderFile_slice.rdbuf()->pubsetbuf(io_buffer_slice.dataPtr(),
+ io_buffer_slice.size());
+ std::string HeaderFileName_slice(WarpX::lab_data_directory+
+ "/slices/Header");
+ HeaderFile_slice.open(HeaderFileName_slice.c_str(),
+ std::ofstream::out |
+ std::ofstream::trunc |
+ std::ofstream::binary);
+
+ if (!HeaderFile_slice.good())
+ FileOpenFailed(HeaderFileName_slice);
+
+ HeaderFile_slice.precision(17);
+
+ HeaderFile_slice << N_slice_snapshots_ << "\n";
+ HeaderFile_slice << dt_slice_snapshots_lab_ << "\n";
+ HeaderFile_slice << gamma_boost_ << "\n";
+ HeaderFile_slice << beta_boost_ << "\n";
+
+ }
+
}
+
+
}
-BoostedFrameDiagnostic::LabSnapShot::
-LabSnapShot(Real t_lab_in, Real t_boost, RealBox prob_domain_lab,
- IntVect prob_ncells_lab,
- int ncomp_to_dump,
- std::vector<std::string> mesh_field_names,
- int file_num_in,
- const BoostedFrameDiagnostic& bfd)
- : t_lab(t_lab_in),
- prob_domain_lab_(prob_domain_lab),
- prob_ncells_lab_(prob_ncells_lab),
- ncomp_to_dump_(ncomp_to_dump),
- mesh_field_names_(mesh_field_names),
- file_num(file_num_in),
- my_bfd(bfd)
+LabFrameSnapShot::
+LabFrameSnapShot(Real t_lab_in, Real t_boost, Real inv_gamma_boost_in,
+ Real inv_beta_boost_in, Real dz_lab_in, RealBox prob_domain_lab,
+ IntVect prob_ncells_lab, int ncomp_to_dump,
+ std::vector<std::string> mesh_field_names,
+ amrex::RealBox diag_domain_lab, Box diag_box, int file_num_in)
{
- Real zmin_lab = prob_domain_lab_.lo(AMREX_SPACEDIM-1);
- current_z_lab = 0.0;
- current_z_boost = 0.0;
- updateCurrentZPositions(t_boost, my_bfd.inv_gamma_boost_, my_bfd.inv_beta_boost_);
- initial_i = (current_z_lab - zmin_lab) / my_bfd.dz_lab_;
- file_name = Concatenate(WarpX::lab_data_directory + "/snapshot", file_num, 5);
+ t_lab = t_lab_in;
+ dz_lab_ = dz_lab_in;
+ inv_gamma_boost_ = inv_gamma_boost_in;
+ inv_beta_boost_ = inv_beta_boost_in;
+ prob_domain_lab_ = prob_domain_lab;
+ prob_ncells_lab_ = prob_ncells_lab;
+ diag_domain_lab_ = diag_domain_lab;
+ buff_box_ = diag_box;
+ ncomp_to_dump_ = ncomp_to_dump;
+ mesh_field_names_ = mesh_field_names;
+ file_num = file_num_in;
+ current_z_lab = 0.0;
+ current_z_boost = 0.0;
+ updateCurrentZPositions(t_boost, inv_gamma_boost_, inv_beta_boost_);
+ Real zmin_lab = prob_domain_lab_.lo(AMREX_SPACEDIM-1);
+ initial_i = (current_z_lab - zmin_lab) / dz_lab_ ;
+ file_name = Concatenate(WarpX::lab_data_directory + "/snapshots/snapshot",
+ file_num, 5);
+ createLabFrameDirectories();
+ buff_counter_ = 0;
+ if (WarpX::do_boosted_frame_fields) data_buffer_.reset(nullptr);
+}
+
+void
+LabFrameDiag::
+updateCurrentZPositions(Real t_boost, Real inv_gamma, Real inv_beta)
+{
+ current_z_boost = (t_lab*inv_gamma - t_boost)*PhysConst::c*inv_beta;
+ current_z_lab = (t_lab - t_boost*inv_gamma)*PhysConst::c*inv_beta;
+}
+void
+LabFrameDiag::
+createLabFrameDirectories() {
#ifdef WARPX_USE_HDF5
if (ParallelDescriptor::IOProcessor())
{
@@ -974,7 +1160,8 @@ LabSnapShot(Real t_lab_in, Real t_boost, RealBox prob_domain_lab,
{
if (WarpX::do_boosted_frame_fields)
{
- for (int comp = 0; comp < ncomp_to_dump; ++comp) {
+ const auto lo = lbound(buff_box_);
+ for (int comp = 0; comp < ncomp_to_dump_; ++comp) {
output_create_field(file_name, mesh_field_names_[comp],
prob_ncells_lab_[0],
#if ( AMREX_SPACEDIM == 3 )
@@ -1036,19 +1223,12 @@ LabSnapShot(Real t_lab_in, Real t_boost, RealBox prob_domain_lab,
#endif
ParallelDescriptor::Barrier();
- writeSnapShotHeader();
-}
-
-void
-BoostedFrameDiagnostic::LabSnapShot::
-updateCurrentZPositions(Real t_boost, Real inv_gamma, Real inv_beta) {
- current_z_boost = (t_lab*inv_gamma - t_boost)*PhysConst::c*inv_beta;
- current_z_lab = (t_lab - t_boost*inv_gamma)*PhysConst::c*inv_beta;
+ writeLabFrameHeader();
}
void
-BoostedFrameDiagnostic::LabSnapShot::
-writeSnapShotHeader() {
+LabFrameDiag::
+writeLabFrameHeader() {
#ifndef WARPX_USE_HDF5
if (ParallelDescriptor::IOProcessor()) {
VisMF::IO_Buffer io_buffer(VisMF::IO_Buffer_Size);
@@ -1072,17 +1252,17 @@ writeSnapShotHeader() {
<< prob_ncells_lab_[AMREX_SPACEDIM-1] <<'\n';
// Write domain physical boundaries
// domain lower bound
- HeaderFile << prob_domain_lab_.lo(0) << ' '
+ HeaderFile << diag_domain_lab_.lo(0) << ' '
#if ( AMREX_SPACEDIM==3 )
- << prob_domain_lab_.lo(1) << ' '
+ << diag_domain_lab_.lo(1) << ' '
#endif
- << prob_domain_lab_.lo(AMREX_SPACEDIM-1) <<'\n';
+ << diag_domain_lab_.lo(AMREX_SPACEDIM-1) <<'\n';
// domain higher bound
- HeaderFile << prob_domain_lab_.hi(0) << ' '
+ HeaderFile << diag_domain_lab_.hi(0) << ' '
#if ( AMREX_SPACEDIM==3 )
- << prob_domain_lab_.hi(1) << ' '
+ << diag_domain_lab_.hi(1) << ' '
#endif
- << prob_domain_lab_.hi(AMREX_SPACEDIM-1) <<'\n';
+ << diag_domain_lab_.hi(AMREX_SPACEDIM-1) <<'\n';
// List of fields dumped to file
for (int i=0; i<ncomp_to_dump_; i++)
{
@@ -1091,4 +1271,197 @@ writeSnapShotHeader() {
HeaderFile << "\n";
}
#endif
+
+}
+
+
+LabFrameSlice::
+LabFrameSlice(Real t_lab_in, Real t_boost, Real inv_gamma_boost_in,
+ Real inv_beta_boost_in, Real dz_lab_in, RealBox prob_domain_lab,
+ IntVect prob_ncells_lab, int ncomp_to_dump,
+ std::vector<std::string> mesh_field_names,
+ RealBox diag_domain_lab, Box diag_box, int file_num_in,
+ amrex::Real cell_dx, amrex::Real cell_dy)
+{
+ t_lab = t_lab_in;
+ prob_domain_lab_ = prob_domain_lab;
+ prob_ncells_lab_ = prob_ncells_lab;
+ diag_domain_lab_ = diag_domain_lab;
+ buff_box_ = diag_box;
+ ncomp_to_dump_ = ncomp_to_dump;
+ mesh_field_names_ = mesh_field_names;
+ file_num = file_num_in;
+ current_z_lab = 0.0;
+ current_z_boost = 0.0;
+ updateCurrentZPositions(t_boost, inv_gamma_boost_, inv_beta_boost_);
+ Real zmin_lab = prob_domain_lab_.lo(AMREX_SPACEDIM-1);
+ initial_i = (current_z_lab - zmin_lab)/dz_lab_;
+ file_name = Concatenate(WarpX::lab_data_directory+"/slices/slice",file_num,5);
+ createLabFrameDirectories();
+ buff_counter_ = 0;
+ dx_ = cell_dx;
+ dy_ = cell_dy;
+
+ if (WarpX::do_boosted_frame_fields) data_buffer_.reset(nullptr);
+}
+
+void
+LabFrameSnapShot::
+AddDataToBuffer( MultiFab& tmp, int k_lab,
+ amrex::Gpu::ManagedDeviceVector<int> map_actual_fields_to_dump)
+{
+ const int ncomp_to_dump = map_actual_fields_to_dump.size();
+ MultiFab& buf = *data_buffer_;
+ for (MFIter mfi(tmp, TilingIfNotGPU()); mfi.isValid(); ++mfi) {
+ Array4<Real> tmp_arr = tmp[mfi].array();
+ Array4<Real> buf_arr = buf[mfi].array();
+ // For 3D runs, tmp is a 2D (x,y) multifab that contains only
+ // slice to write to file
+ const Box& bx = mfi.tilebox();
+ const auto field_map_ptr = map_actual_fields_to_dump.dataPtr();
+ ParallelFor(bx, ncomp_to_dump,
+ [=] AMREX_GPU_DEVICE(int i, int j, int k, int n)
+ {
+ const int icomp = field_map_ptr[n];
+#if (AMREX_SPACEDIM == 3)
+ buf_arr(i,j,k_lab,n) = tmp_arr(i,j,k,icomp);
+#else
+ buf_arr(i,k_lab,k,n) = tmp_arr(i,j,k,icomp);
+#endif
+ }
+ );
+ }
+}
+
+
+void
+LabFrameSlice::
+AddDataToBuffer( MultiFab& tmp, int k_lab,
+ amrex::Gpu::ManagedDeviceVector<int> map_actual_fields_to_dump)
+{
+ const int ncomp_to_dump = map_actual_fields_to_dump.size();
+ MultiFab& buf = *data_buffer_;
+ for (MFIter mfi(tmp, TilingIfNotGPU()); mfi.isValid(); ++mfi)
+ {
+ Box& bx = buff_box_;
+ const Box& bx_bf = mfi.tilebox();
+ bx.setSmall(AMREX_SPACEDIM-1,bx_bf.smallEnd(AMREX_SPACEDIM-1));
+ bx.setBig(AMREX_SPACEDIM-1,bx_bf.bigEnd(AMREX_SPACEDIM-1));
+ Array4<Real> tmp_arr = tmp[mfi].array();
+ Array4<Real> buf_arr = buf[mfi].array();
+ const auto field_map_ptr = map_actual_fields_to_dump.dataPtr();
+ ParallelFor(bx, ncomp_to_dump,
+ [=] AMREX_GPU_DEVICE(int i, int j, int k, int n)
+ {
+ const int icomp = field_map_ptr[n];
+#if (AMREX_SPACEDIM == 3)
+ buf_arr(i,j,k_lab,n) = tmp_arr(i,j,k,icomp);
+#else
+ buf_arr(i,k_lab,k,n) = tmp_arr(i,j,k,icomp);
+#endif
+ });
+ }
+
+}
+
+
+void
+LabFrameSnapShot::
+AddPartDataToParticleBuffer(
+ Vector<WarpXParticleContainer::DiagnosticParticleData> tmp_particle_buffer,
+ int nspeciesBoostedFrame) {
+ for (int isp = 0; isp < nspeciesBoostedFrame; ++isp) {
+ auto np = tmp_particle_buffer[isp].GetRealData(DiagIdx::w).size();
+ if (np == 0) return;
+
+ particles_buffer_[isp].GetRealData(DiagIdx::w).insert(
+ particles_buffer_[isp].GetRealData(DiagIdx::w).end(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::w).begin(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::w).end());
+
+ particles_buffer_[isp].GetRealData(DiagIdx::x).insert(
+ particles_buffer_[isp].GetRealData(DiagIdx::x).end(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::x).begin(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::x).end());
+
+ particles_buffer_[isp].GetRealData(DiagIdx::y).insert(
+ particles_buffer_[isp].GetRealData(DiagIdx::y).end(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::y).begin(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::y).end());
+
+ particles_buffer_[isp].GetRealData(DiagIdx::z).insert(
+ particles_buffer_[isp].GetRealData(DiagIdx::z).end(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::z).begin(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::z).end());
+
+ particles_buffer_[isp].GetRealData(DiagIdx::ux).insert(
+ particles_buffer_[isp].GetRealData(DiagIdx::ux).end(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::ux).begin(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::ux).end());
+
+ particles_buffer_[isp].GetRealData(DiagIdx::uy).insert(
+ particles_buffer_[isp].GetRealData(DiagIdx::uy).end(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::uy).begin(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::uy).end());
+
+ particles_buffer_[isp].GetRealData(DiagIdx::uz).insert(
+ particles_buffer_[isp].GetRealData(DiagIdx::uz).end(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::uz).begin(),
+ tmp_particle_buffer[isp].GetRealData(DiagIdx::uz).end());
+ }
+}
+
+void
+LabFrameSlice::
+AddPartDataToParticleBuffer(
+ Vector<WarpXParticleContainer::DiagnosticParticleData> tmp_particle_buffer,
+ int nSpeciesBoostedFrame) {
+ for (int isp = 0; isp < nSpeciesBoostedFrame; ++isp) {
+ auto np = tmp_particle_buffer[isp].GetRealData(DiagIdx::w).size();
+
+ if (np == 0) return;
+
+ auto const& wpc = tmp_particle_buffer[isp].GetRealData(DiagIdx::w);
+ auto const& xpc = tmp_particle_buffer[isp].GetRealData(DiagIdx::x);
+ auto const& ypc = tmp_particle_buffer[isp].GetRealData(DiagIdx::y);
+ auto const& zpc = tmp_particle_buffer[isp].GetRealData(DiagIdx::z);
+ auto const& uxpc = tmp_particle_buffer[isp].GetRealData(DiagIdx::ux);
+ auto const& uypc = tmp_particle_buffer[isp].GetRealData(DiagIdx::uy);
+ auto const& uzpc = tmp_particle_buffer[isp].GetRealData(DiagIdx::uz);
+
+ particles_buffer_[isp].resize(np);
+ auto wpc_buff = particles_buffer_[isp].GetRealData(DiagIdx::w);
+ auto xpc_buff = particles_buffer_[isp].GetRealData(DiagIdx::x);
+ auto ypc_buff = particles_buffer_[isp].GetRealData(DiagIdx::y);
+ auto zpc_buff = particles_buffer_[isp].GetRealData(DiagIdx::z);
+ auto uxpc_buff = particles_buffer_[isp].GetRealData(DiagIdx::ux);
+ auto uypc_buff = particles_buffer_[isp].GetRealData(DiagIdx::uy);
+ auto uzpc_buff = particles_buffer_[isp].GetRealData(DiagIdx::uz);
+
+
+ int partcounter = 0;
+ for (int i = 0; i < np; ++i)
+ {
+ if( xpc[i] >= (diag_domain_lab_.lo(0)-dx_) &&
+ xpc[i] <= (diag_domain_lab_.hi(0)+dx_) ) {
+ #if (AMREX_SPACEDIM == 3)
+ if( ypc[i] >= (diag_domain_lab_.lo(1)-dy_) &&
+ ypc[i] <= (diag_domain_lab_.hi(1) + dy_))
+ #endif
+ {
+ wpc_buff[partcounter] = wpc[i];
+ xpc_buff[partcounter] = xpc[i];
+ ypc_buff[partcounter] = ypc[i];
+ zpc_buff[partcounter] = zpc[i];
+ uxpc_buff[partcounter] = uxpc[i];
+ uypc_buff[partcounter] = uypc[i];
+ uzpc_buff[partcounter] = uzpc[i];
+ ++partcounter;
+ }
+ }
+ }
+
+ particles_buffer_[isp].resize(partcounter);
+
+ }
}