aboutsummaryrefslogtreecommitdiff
path: root/Source
diff options
context:
space:
mode:
Diffstat (limited to 'Source')
-rw-r--r--Source/CustomDensityProb.cpp55
-rw-r--r--Source/ParticleContainer.H2
-rw-r--r--Source/ParticleContainer.cpp8
-rw-r--r--Source/PhysicalParticleContainer.cpp4
-rw-r--r--Source/WarpX.H12
-rw-r--r--Source/WarpXBoostedFrameDiagnostic.H4
-rw-r--r--Source/WarpXBoostedFrameDiagnostic.cpp30
-rw-r--r--Source/WarpXEvolve.cpp7
-rw-r--r--Source/WarpXIO.cpp2
-rw-r--r--Source/WarpXInitData.cpp8
-rw-r--r--Source/WarpXWrappers.cpp72
11 files changed, 175 insertions, 29 deletions
diff --git a/Source/CustomDensityProb.cpp b/Source/CustomDensityProb.cpp
index 1df3d75ad..2f6005bc2 100644
--- a/Source/CustomDensityProb.cpp
+++ b/Source/CustomDensityProb.cpp
@@ -1,13 +1,52 @@
#include <PlasmaInjector.H>
-#include <iostream>
+using namespace amrex;
-#include <AMReX.H>
+///
+/// This "custom" density profile just does constant
+///
+Real CustomDensityProfile::getDensity(Real x, Real y, Real z) const {
+ const Real on_axis_density = params[0];
+ const Real plasma_zmin = params[1];
+ const Real plasma_zmax = params[2];
+ const Real plasma_lramp_start = params[3];
+ const Real plasma_lramp_end = params[4];
+ const Real plasma_rcap = params[5];
+ const Real plasma_rdownramp = params[6];
+ const Real plasma_rchannel = params[7];
+ static const Real re = 2.8178403227e-15; // Electron classical radius
+ static const Real pi = 3.14159265359;
-amrex::Real CustomDensityProfile::getDensity(amrex::Real x,
- amrex::Real y,
- amrex::Real z) const
-{
- amrex::Abort("If running with a custom density profile, you must supply a CustomDensityProb.cpp file");
- return 0.0;
+ Real r2 = x*x + y*y;
+ Real r = std::sqrt( r2 );
+
+ // Transverse part of the profile
+ Real nr;
+ if (r<plasma_rcap) {
+ nr = 1. + 1./(pi*on_axis_density*re) * r2/pow(plasma_rchannel, 4);
+ } else {
+ nr = 1. + 1./(pi*on_axis_density*re) *
+ pow(plasma_rcap, 2)/pow(plasma_rchannel, 4) *
+ (plasma_rcap+plasma_rdownramp-r)/plasma_rdownramp;
+ }
+ // Longitudinal part of the profile
+ Real nz;
+ if (z<plasma_zmin) {
+ nz = 0;
+ } else if (z<plasma_zmin+plasma_lramp_start) {
+ nz = (z-plasma_zmin)/plasma_lramp_start;
+ } else if (z<plasma_zmax-plasma_lramp_end) {
+ nz = 1.;
+ } else if (z<plasma_zmax){
+ nz = -(z-plasma_zmax)/plasma_lramp_end;
+ } else {
+ nz = 0;
+ }
+ // Combine and saturate profile
+ Real n = nr*nz;
+ if (n > 4.) {
+ n = 4.;
+ }
+
+ return on_axis_density*n;
}
diff --git a/Source/ParticleContainer.H b/Source/ParticleContainer.H
index a1402d4df..9767edd59 100644
--- a/Source/ParticleContainer.H
+++ b/Source/ParticleContainer.H
@@ -120,6 +120,8 @@ public:
void Redistribute ();
+ void RedistributeLocal ();
+
amrex::Vector<long> NumberOfParticlesInGrid(int lev) const;
void Increment (amrex::MultiFab& mf, int lev);
diff --git a/Source/ParticleContainer.cpp b/Source/ParticleContainer.cpp
index 55889cf8d..168b9167c 100644
--- a/Source/ParticleContainer.cpp
+++ b/Source/ParticleContainer.cpp
@@ -186,6 +186,14 @@ MultiParticleContainer::Redistribute ()
}
}
+void
+MultiParticleContainer::RedistributeLocal ()
+{
+ for (auto& pc : allcontainers) {
+ pc->Redistribute(0, 0, 0, true);
+ }
+}
+
Vector<long>
MultiParticleContainer::NumberOfParticlesInGrid(int lev) const
{
diff --git a/Source/PhysicalParticleContainer.cpp b/Source/PhysicalParticleContainer.cpp
index e5f6b0c82..050e50daa 100644
--- a/Source/PhysicalParticleContainer.cpp
+++ b/Source/PhysicalParticleContainer.cpp
@@ -225,9 +225,9 @@ PhysicalParticleContainer::AddPlasma(int lev, RealBox part_realbox )
// and that the boost is along z)
Real t = WarpX::GetInstance().gett_new(lev);
Real v_boost = WarpX::beta_boost*PhysConst::c;
- Real z_lab = WarpX::gamma_boost*( z - v_boost*t );
+ Real z_lab = WarpX::gamma_boost*( z + v_boost*t );
plasma_injector->getMomentum(u, x, y, z_lab);
- dens = plasma_injector->getDensity(x, y, z);
+ dens = plasma_injector->getDensity(x, y, z_lab);
// Perform Lorentz transform
// (Assumes that the plasma has a low velocity)
u[2] = WarpX::gamma_boost * ( u[2] - v_boost );
diff --git a/Source/WarpX.H b/Source/WarpX.H
index d0876afd2..35cf93cc7 100644
--- a/Source/WarpX.H
+++ b/Source/WarpX.H
@@ -90,8 +90,16 @@ public:
static amrex::Vector<amrex::Real> boost_direction;
const amrex::MultiFab& getcurrent (int lev, int direction) {return *current_fp[lev][direction];}
- const amrex::MultiFab& getEfield (int lev, int direction) {return *Efield_fp[lev][direction];}
- const amrex::MultiFab& getBfield (int lev, int direction) {return *Bfield_fp[lev][direction];}
+ const amrex::MultiFab& getEfield (int lev, int direction) {return *Efield_aux[lev][direction];}
+ const amrex::MultiFab& getBfield (int lev, int direction) {return *Bfield_aux[lev][direction];}
+
+ const amrex::MultiFab& getcurrent_cp (int lev, int direction) {return *current_cp[lev][direction];}
+ const amrex::MultiFab& getEfield_cp (int lev, int direction) {return *Efield_cp[lev][direction];}
+ const amrex::MultiFab& getBfield_cp (int lev, int direction) {return *Bfield_cp[lev][direction];}
+
+ const amrex::MultiFab& getcurrent_fp (int lev, int direction) {return *current_fp[lev][direction];}
+ const amrex::MultiFab& getEfield_fp (int lev, int direction) {return *Efield_fp[lev][direction];}
+ const amrex::MultiFab& getBfield_fp (int lev, int direction) {return *Bfield_fp[lev][direction];}
static amrex::MultiFab* getCosts (int lev) {
if (m_instance) {
diff --git a/Source/WarpXBoostedFrameDiagnostic.H b/Source/WarpXBoostedFrameDiagnostic.H
index a5f0bd7ff..96a77f182 100644
--- a/Source/WarpXBoostedFrameDiagnostic.H
+++ b/Source/WarpXBoostedFrameDiagnostic.H
@@ -67,8 +67,8 @@ class BoostedFrameDiagnostic {
int boost_direction_;
amrex::Vector<std::unique_ptr<amrex::MultiFab> > data_buffer_;
- int num_buffer_ = 32;
- int max_box_size_ = 64;
+ int num_buffer_ = 256;
+ int max_box_size_ = 256;
amrex::Vector<int> buff_counter_;
amrex::Vector<LabSnapShot> snapshots_;
diff --git a/Source/WarpXBoostedFrameDiagnostic.cpp b/Source/WarpXBoostedFrameDiagnostic.cpp
index d58ed0be7..a617961ee 100644
--- a/Source/WarpXBoostedFrameDiagnostic.cpp
+++ b/Source/WarpXBoostedFrameDiagnostic.cpp
@@ -46,7 +46,10 @@ BoostedFrameDiagnostic::
writeLabFrameData(const MultiFab& cell_centered_data, const Geometry& geom, Real t_boost) {
BL_PROFILE("BoostedFrameDiagnostic::writeLabFrameData");
-
+
+ VisMF::Header::Version current_version = VisMF::GetHeaderVersion();
+ VisMF::SetHeaderVersion(amrex::VisMF::Header::NoFabHeader_v1);
+
const RealBox& domain_z_boost = geom.ProbDomain();
const Real zlo_boost = domain_z_boost.lo(boost_direction_);
const Real zhi_boost = domain_z_boost.hi(boost_direction_);
@@ -119,6 +122,8 @@ writeLabFrameData(const MultiFab& cell_centered_data, const Geometry& geom, Real
buff_counter_[i] = 0;
}
}
+
+ VisMF::SetHeaderVersion(current_version);
}
void
@@ -167,20 +172,23 @@ LabSnapShot(Real t_lab_in, Real zmin_lab_in,
current_z_boost = 0.0;
file_name = Concatenate("lab_frame_data/snapshot", file_num, 5);
- const int nlevels = 1;
- const std::string level_prefix = "Level_";
-
- if (!UtilCreateDirectory(file_name, 0755))
- CreateDirectoryFailed(file_name);
- for(int i(0); i < nlevels; ++i) {
- const std::string &fullpath = LevelFullPath(i, file_name);
- if (!UtilCreateDirectory(fullpath, 0755))
- CreateDirectoryFailed(fullpath);
+ if (ParallelDescriptor::IOProcessor()) {
+
+ const int nlevels = 1;
+ const std::string level_prefix = "Level_";
+
+ if (!UtilCreateDirectory(file_name, 0755))
+ CreateDirectoryFailed(file_name);
+ for(int i(0); i < nlevels; ++i) {
+ const std::string &fullpath = LevelFullPath(i, file_name);
+ if (!UtilCreateDirectory(fullpath, 0755))
+ CreateDirectoryFailed(fullpath);
+ }
}
ParallelDescriptor::Barrier();
- writeSnapShotHeader();
+ if (ParallelDescriptor::IOProcessor()) writeSnapShotHeader();
}
void
diff --git a/Source/WarpXEvolve.cpp b/Source/WarpXEvolve.cpp
index 6ce7db0ef..cda70dd53 100644
--- a/Source/WarpXEvolve.cpp
+++ b/Source/WarpXEvolve.cpp
@@ -260,7 +260,12 @@ WarpX::EvolveEM (int numsteps)
// We might need to move j because we are going to make a plotfile.
MoveWindow(move_j);
- mypc->Redistribute(); // Redistribute particles
+ if (max_level == 0) {
+ mypc->RedistributeLocal();
+ }
+ else {
+ mypc->Redistribute();
+ }
amrex::Print()<< "STEP " << step+1 << " ends." << " TIME = " << cur_time
<< " DT = " << dt[0] << "\n";
diff --git a/Source/WarpXIO.cpp b/Source/WarpXIO.cpp
index b289fdde9..212d69a14 100644
--- a/Source/WarpXIO.cpp
+++ b/Source/WarpXIO.cpp
@@ -404,7 +404,7 @@ WarpX::GetCellCenteredData() {
DistributionMap(lev),
nc, ng) );
- Array<const MultiFab*> srcmf(BL_SPACEDIM);
+ Vector<const MultiFab*> srcmf(BL_SPACEDIM);
int dcomp = 0;
// first the electric field
diff --git a/Source/WarpXInitData.cpp b/Source/WarpXInitData.cpp
index df01afe88..0d6c35a4d 100644
--- a/Source/WarpXInitData.cpp
+++ b/Source/WarpXInitData.cpp
@@ -54,8 +54,12 @@ WarpX::InitDiagnostics () {
const Real* current_hi = geom[0].ProbHi();
Real dt_boost = dt[0];
- myBFD.reset(new BoostedFrameDiagnostic(current_lo[moving_window_dir],
- current_hi[moving_window_dir],
+ // Find the positions of the lab-frame box that corresponds to the boosted-frame box at t=0
+ Real zmin_lab = current_lo[moving_window_dir]/( (1.+beta_boost)*gamma_boost );
+ Real zmax_lab = current_hi[moving_window_dir]/( (1.+beta_boost)*gamma_boost );
+
+ myBFD.reset(new BoostedFrameDiagnostic(zmin_lab,
+ zmax_lab,
moving_window_v, dt_snapshots_lab,
num_snapshots_lab, gamma_boost, dt_boost,
moving_window_dir));
diff --git a/Source/WarpXWrappers.cpp b/Source/WarpXWrappers.cpp
index 54aeedf35..d106cfca7 100644
--- a/Source/WarpXWrappers.cpp
+++ b/Source/WarpXWrappers.cpp
@@ -192,6 +192,30 @@ extern "C"
return getMultiFabLoVects(mf, return_size, ngrow);
}
+ double** warpx_getEfieldCP(int lev, int direction,
+ int *return_size, int *ngrow, int **shapes) {
+ auto & mf = WarpX::GetInstance().getEfield_cp(lev, direction);
+ return getMultiFabPointers(mf, return_size, ngrow, shapes);
+ }
+
+ int* warpx_getEfieldCPLoVects(int lev, int direction,
+ int *return_size, int *ngrow) {
+ auto & mf = WarpX::GetInstance().getEfield_cp(lev, direction);
+ return getMultiFabLoVects(mf, return_size, ngrow);
+ }
+
+ double** warpx_getEfieldFP(int lev, int direction,
+ int *return_size, int *ngrow, int **shapes) {
+ auto & mf = WarpX::GetInstance().getEfield_fp(lev, direction);
+ return getMultiFabPointers(mf, return_size, ngrow, shapes);
+ }
+
+ int* warpx_getEfieldFPLoVects(int lev, int direction,
+ int *return_size, int *ngrow) {
+ auto & mf = WarpX::GetInstance().getEfield_fp(lev, direction);
+ return getMultiFabLoVects(mf, return_size, ngrow);
+ }
+
double** warpx_getBfield(int lev, int direction,
int *return_size, int *ngrow, int **shapes) {
auto & mf = WarpX::GetInstance().getBfield(lev, direction);
@@ -204,6 +228,30 @@ extern "C"
return getMultiFabLoVects(mf, return_size, ngrow);
}
+ double** warpx_getBfieldCP(int lev, int direction,
+ int *return_size, int *ngrow, int **shapes) {
+ auto & mf = WarpX::GetInstance().getBfield_cp(lev, direction);
+ return getMultiFabPointers(mf, return_size, ngrow, shapes);
+ }
+
+ int* warpx_getBfieldCPLoVects(int lev, int direction,
+ int *return_size, int *ngrow) {
+ auto & mf = WarpX::GetInstance().getBfield_cp(lev, direction);
+ return getMultiFabLoVects(mf, return_size, ngrow);
+ }
+
+ double** warpx_getBfieldFP(int lev, int direction,
+ int *return_size, int *ngrow, int **shapes) {
+ auto & mf = WarpX::GetInstance().getBfield_fp(lev, direction);
+ return getMultiFabPointers(mf, return_size, ngrow, shapes);
+ }
+
+ int* warpx_getBfieldFPLoVects(int lev, int direction,
+ int *return_size, int *ngrow) {
+ auto & mf = WarpX::GetInstance().getBfield_fp(lev, direction);
+ return getMultiFabLoVects(mf, return_size, ngrow);
+ }
+
double** warpx_getCurrentDensity(int lev, int direction,
int *return_size, int *ngrow, int **shapes) {
auto & mf = WarpX::GetInstance().getcurrent(lev, direction);
@@ -216,6 +264,30 @@ extern "C"
return getMultiFabLoVects(mf, return_size, ngrow);
}
+ double** warpx_getCurrentDensityCP(int lev, int direction,
+ int *return_size, int *ngrow, int **shapes) {
+ auto & mf = WarpX::GetInstance().getcurrent_cp(lev, direction);
+ return getMultiFabPointers(mf, return_size, ngrow, shapes);
+ }
+
+ int* warpx_getCurrentDensityCPLoVects(int lev, int direction,
+ int *return_size, int *ngrow) {
+ auto & mf = WarpX::GetInstance().getcurrent_cp(lev, direction);
+ return getMultiFabLoVects(mf, return_size, ngrow);
+ }
+
+ double** warpx_getCurrentDensityFP(int lev, int direction,
+ int *return_size, int *ngrow, int **shapes) {
+ auto & mf = WarpX::GetInstance().getcurrent_fp(lev, direction);
+ return getMultiFabPointers(mf, return_size, ngrow, shapes);
+ }
+
+ int* warpx_getCurrentDensityFPLoVects(int lev, int direction,
+ int *return_size, int *ngrow) {
+ auto & mf = WarpX::GetInstance().getcurrent_fp(lev, direction);
+ return getMultiFabLoVects(mf, return_size, ngrow);
+ }
+
double** warpx_getParticleStructs(int speciesnumber,
int* num_tiles, int** particles_per_tile) {
auto & mypc = WarpX::GetInstance().GetPartContainer();
generated by cgit v1.2.3 (git 2.46.0) at 2025-08-11 23:33:23 +0000