aboutsummaryrefslogtreecommitdiff
path: root/Source/Initialization
diff options
context:
space:
mode:
Diffstat (limited to 'Source/Initialization')
-rw-r--r--Source/Initialization/CustomDensityProb.H2
-rw-r--r--Source/Initialization/InjectorDensity.H17
-rw-r--r--Source/Initialization/InjectorDensity.cpp3
-rw-r--r--Source/Initialization/InjectorMomentum.H15
-rw-r--r--Source/Initialization/InjectorMomentum.cpp1
-rw-r--r--Source/Initialization/InjectorPosition.H14
-rw-r--r--Source/Initialization/PlasmaInjector.H16
-rw-r--r--Source/Initialization/PlasmaInjector.cpp17
-rw-r--r--Source/Initialization/WarpXInitData.cpp51
9 files changed, 76 insertions, 60 deletions
diff --git a/Source/Initialization/CustomDensityProb.H b/Source/Initialization/CustomDensityProb.H
index b00830e6c..c5c159ee8 100644
--- a/Source/Initialization/CustomDensityProb.H
+++ b/Source/Initialization/CustomDensityProb.H
@@ -27,7 +27,7 @@ struct InjectorDensityCustom
}
}
- // Return density at given position, using user-defined parameters
+ // Return density at given position, using user-defined parameters
// stored in p.
AMREX_GPU_HOST_DEVICE
amrex::Real
diff --git a/Source/Initialization/InjectorDensity.H b/Source/Initialization/InjectorDensity.H
index b7f5c26eb..7e61ae27d 100644
--- a/Source/Initialization/InjectorDensity.H
+++ b/Source/Initialization/InjectorDensity.H
@@ -1,12 +1,13 @@
#ifndef INJECTOR_DENSITY_H_
#define INJECTOR_DENSITY_H_
-#include <AMReX_Gpu.H>
-#include <AMReX_Dim3.H>
#include <GpuParser.H>
#include <CustomDensityProb.H>
#include <WarpXConst.H>
+#include <AMReX_Gpu.H>
+#include <AMReX_Dim3.H>
+
// struct whose getDensity returns constant density.
struct InjectorDensityConstant
{
@@ -67,18 +68,20 @@ struct InjectorDensityPredefined
amrex::Real kp = PhysConst::q_e/PhysConst::c
*std::sqrt( n0/(PhysConst::m_e*PhysConst::ep0) );
+ // Longitudinal profile, normalized to 1
if ((z-z_start)>=0 and
(z-z_start)<ramp_up ) {
- n = (z-z_start)/ramp_up;
+ n = 0.5*(1.-std::cos(MathConst::pi*(z-z_start)/ramp_up));
} else if ((z-z_start)>=ramp_up and
(z-z_start)< ramp_up+plateau ) {
n = 1.;
} else if ((z-z_start)>=ramp_up+plateau and
(z-z_start)< ramp_up+plateau+ramp_down) {
- n = 1.-((z-z_start)-ramp_up-plateau)/ramp_down;
+ n = 0.5*(1.+std::cos(MathConst::pi*((z-z_start)-ramp_up-plateau)/ramp_down));
} else {
n = 0.;
}
+ // Multiply by transverse profile, and physical density
n *= n0*(1.+4.*(x*x+y*y)/(kp*kp*rc*rc*rc*rc));
return n;
}
@@ -94,9 +97,9 @@ private:
amrex::Real* p;
};
-// Base struct for density injector.
-// InjectorDensity contains a union (called Object) that holds any one
-// instance of:
+// Base struct for density injector.
+// InjectorDensity contains a union (called Object) that holds any one
+// instance of:
// - InjectorDensityConstant : to generate constant density;
// - InjectorDensityParser : to generate density from parser;
// - InjectorDensityCustom : to generate density from custom profile;
diff --git a/Source/Initialization/InjectorDensity.cpp b/Source/Initialization/InjectorDensity.cpp
index 54df4b14d..ddf6d1a18 100644
--- a/Source/Initialization/InjectorDensity.cpp
+++ b/Source/Initialization/InjectorDensity.cpp
@@ -1,3 +1,4 @@
+#include <InjectorDensity.H>
#include <PlasmaInjector.H>
using namespace amrex;
@@ -48,7 +49,7 @@ InjectorDensityPredefined::InjectorDensityPredefined (
ParmParse pp(a_species_name);
std::vector<amrex::Real> v;
- // Read parameters for the predefined plasma profile,
+ // Read parameters for the predefined plasma profile,
// and store them in managed memory
pp.getarr("predefined_profile_params", v);
p = static_cast<amrex::Real*>
diff --git a/Source/Initialization/InjectorMomentum.H b/Source/Initialization/InjectorMomentum.H
index 399ee7759..ba0d26fc4 100644
--- a/Source/Initialization/InjectorMomentum.H
+++ b/Source/Initialization/InjectorMomentum.H
@@ -1,10 +1,11 @@
#ifndef INJECTOR_MOMENTUM_H_
#define INJECTOR_MOMENTUM_H_
+#include <CustomMomentumProb.H>
+#include <GpuParser.H>
+
#include <AMReX_Gpu.H>
#include <AMReX_Dim3.H>
-#include <GpuParser.H>
-#include <CustomMomentumProb.H>
// struct whose getMomentum returns constant momentum.
struct InjectorMomentumConstant
@@ -22,7 +23,7 @@ private:
amrex::Real m_ux, m_uy, m_uz;
};
-// struct whose getMomentum returns momentum for 1 particle, from random
+// struct whose getMomentum returns momentum for 1 particle, from random
// gaussian distribution.
struct InjectorMomentumGaussian
{
@@ -84,9 +85,9 @@ struct InjectorMomentumParser
GpuParser m_ux_parser, m_uy_parser, m_uz_parser;
};
-// Base struct for momentum injector.
-// InjectorMomentum contains a union (called Object) that holds any one
-// instance of:
+// Base struct for momentum injector.
+// InjectorMomentum contains a union (called Object) that holds any one
+// instance of:
// - InjectorMomentumConstant : to generate constant density;
// - InjectorMomentumGaussian : to generate gaussian distribution;
// - InjectorMomentumRadialExpansion: to generate radial expansion;
@@ -189,7 +190,7 @@ private:
Type type;
// An instance of union Object constructs and stores any one of
- // the objects declared (constant or custom or gaussian or
+ // the objects declared (constant or custom or gaussian or
// radial_expansion or parser).
union Object {
Object (InjectorMomentumConstant*,
diff --git a/Source/Initialization/InjectorMomentum.cpp b/Source/Initialization/InjectorMomentum.cpp
index a197b5bef..8fadf0c4b 100644
--- a/Source/Initialization/InjectorMomentum.cpp
+++ b/Source/Initialization/InjectorMomentum.cpp
@@ -1,3 +1,4 @@
+#include <InjectorMomentum.H>
#include <PlasmaInjector.H>
using namespace amrex;
diff --git a/Source/Initialization/InjectorPosition.H b/Source/Initialization/InjectorPosition.H
index 19bb092dd..6ecae93e0 100644
--- a/Source/Initialization/InjectorPosition.H
+++ b/Source/Initialization/InjectorPosition.H
@@ -25,7 +25,7 @@ struct InjectorPositionRegular
// i_part: particle number within the cell, required to evenly space
// particles within the cell.
- // ref_fac: the number of particles evenly-spaced within a cell
+ // ref_fac: the number of particles evenly-spaced within a cell
// is a_ppc*(ref_fac**AMREX_SPACEDIM).
AMREX_GPU_HOST_DEVICE
amrex::XDim3
@@ -33,7 +33,7 @@ struct InjectorPositionRegular
{
int nx = ref_fac*ppc.x;
int ny = ref_fac*ppc.y;
-#if (AMREX_SPACEDIM == 3)
+#if (defined WARPX_DIM_3D) || (defined WARPX_DIM_RZ)
int nz = ref_fac*ppc.z;
#else
int nz = 1;
@@ -41,15 +41,17 @@ struct InjectorPositionRegular
int ix_part = i_part/(ny*nz); // written this way backward compatibility
int iz_part = (i_part-ix_part*(ny*nz)) / ny;
int iy_part = (i_part-ix_part*(ny*nz)) - ny*iz_part;
- return amrex::XDim3{(0.5+ix_part)/nx, (0.5+iy_part)/ny, (0.5+iz_part) / nz};
+ return amrex::XDim3{(amrex::Real(0.5)+ix_part)/nx,
+ (amrex::Real(0.5)+iy_part)/ny,
+ (amrex::Real(0.5)+iz_part) / nz};
}
private:
amrex::Dim3 ppc;
};
-// Base struct for position injector.
-// InjectorPosition contains a union (called Object) that holds any one
-// instance of:
+// Base struct for position injector.
+// InjectorPosition contains a union (called Object) that holds any one
+// instance of:
// - InjectorPositionRandom : to generate random distribution;
// - InjectorPositionRegular: to generate regular distribution.
// The choice is made at runtime, depending in the constructor called.
diff --git a/Source/Initialization/PlasmaInjector.H b/Source/Initialization/PlasmaInjector.H
index f7e86bff5..56b32c827 100644
--- a/Source/Initialization/PlasmaInjector.H
+++ b/Source/Initialization/PlasmaInjector.H
@@ -5,13 +5,15 @@
#include <InjectorDensity.H>
#include <InjectorMomentum.H>
-#include <array>
-#include <AMReX_Vector.H>
#include <WarpXConst.H>
#include <WarpXParser.H>
+
+#include <AMReX_Vector.H>
#include <AMReX_ParmParse.H>
#include <AMReX_Utility.H>
+#include <array>
+
///
/// The PlasmaInjector class parses and stores information about the plasma
/// type used in the particle container. This information is used to create the
@@ -41,9 +43,9 @@ public:
bool doInjection () const noexcept { return inj_pos != NULL;}
bool add_single_particle = false;
- amrex::Vector<amrex::Real> single_particle_pos;
- amrex::Vector<amrex::Real> single_particle_vel;
- amrex::Real single_particle_weight;
+ amrex::Vector<amrex::ParticleReal> single_particle_pos;
+ amrex::Vector<amrex::ParticleReal> single_particle_vel;
+ amrex::ParticleReal single_particle_weight;
bool gaussian_beam = false;
amrex::Real x_m;
@@ -94,9 +96,9 @@ protected:
std::unique_ptr<InjectorPosition> inj_pos;
std::unique_ptr<InjectorDensity > inj_rho;
std::unique_ptr<InjectorMomentum> inj_mom;
-
+
void parseDensity (amrex::ParmParse& pp);
- void parseMomentum (amrex::ParmParse& pp);
+ void parseMomentum (amrex::ParmParse& pp);
};
#endif
diff --git a/Source/Initialization/PlasmaInjector.cpp b/Source/Initialization/PlasmaInjector.cpp
index 541999789..af04c6b31 100644
--- a/Source/Initialization/PlasmaInjector.cpp
+++ b/Source/Initialization/PlasmaInjector.cpp
@@ -1,13 +1,14 @@
#include "PlasmaInjector.H"
-#include <sstream>
-#include <functional>
-
#include <WarpXConst.H>
#include <WarpX_f.H>
-#include <AMReX.H>
#include <WarpX.H>
+#include <AMReX.H>
+
+#include <sstream>
+#include <functional>
+
using namespace amrex;
namespace {
@@ -45,7 +46,7 @@ namespace {
} else if (name == "m_p"){
return PhysConst::m_p;
} else if (name == "inf"){
- return std::numeric_limits<double>::infinity();
+ return std::numeric_limits<double>::infinity();
} else if (pp.query("mass", result)) {
return result;
} else {
@@ -143,9 +144,11 @@ PlasmaInjector::PlasmaInjector (int ispecies, const std::string& name)
parseDensity(pp);
parseMomentum(pp);
} else if (part_pos_s == "nuniformpercell") {
- num_particles_per_cell_each_dim.resize(3);
+ // Note that for RZ, three numbers are expected, r, theta, and z.
+ // For 2D, only two are expected. The third is overwritten with 1.
+ num_particles_per_cell_each_dim.assign(3, 1);
pp.getarr("num_particles_per_cell_each_dim", num_particles_per_cell_each_dim);
-#if ( AMREX_SPACEDIM == 2 )
+#if WARPX_DIM_XZ
num_particles_per_cell_each_dim[2] = 1;
#endif
// Construct InjectorPosition from InjectorPositionRegular.
diff --git a/Source/Initialization/WarpXInitData.cpp b/Source/Initialization/WarpXInitData.cpp
index 590c11b84..385993f78 100644
--- a/Source/Initialization/WarpXInitData.cpp
+++ b/Source/Initialization/WarpXInitData.cpp
@@ -1,12 +1,12 @@
-#include <AMReX_ParallelDescriptor.H>
-#include <AMReX_ParmParse.H>
-
#include <WarpX.H>
#include <WarpX_f.H>
#include <BilinearFilter.H>
#include <NCIGodfreyFilter.H>
+#include <AMReX_ParallelDescriptor.H>
+#include <AMReX_ParmParse.H>
+
#ifdef BL_USE_SENSEI_INSITU
#include <AMReX_AmrMeshInSituBridge.H>
#endif
@@ -25,11 +25,11 @@ WarpX::InitData ()
}
else
{
- InitFromCheckpoint();
+ InitFromCheckpoint();
if (is_synchronized) {
ComputeDt();
}
- PostRestart();
+ PostRestart();
}
ComputePMLFactors();
@@ -87,12 +87,12 @@ WarpX::InitDiagnostics () {
const Real* current_hi = geom[0].ProbHi();
Real dt_boost = dt[0];
- // 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 );
+ // 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,
+ zmax_lab,
moving_window_v, dt_snapshots_lab,
num_snapshots_lab, gamma_boost,
t_new[0], dt_boost,
@@ -142,6 +142,7 @@ WarpX::InitPML ()
dt[0], nox_fft, noy_fft, noz_fft, do_nodal,
#endif
do_dive_cleaning, do_moving_window,
+ pml_has_particles, do_pml_in_domain,
do_pml_Lo_corrected, do_pml_Hi));
for (int lev = 1; lev <= finest_level; ++lev)
{
@@ -159,6 +160,7 @@ WarpX::InitPML ()
dt[lev], nox_fft, noy_fft, noz_fft, do_nodal,
#endif
do_dive_cleaning, do_moving_window,
+ pml_has_particles, do_pml_in_domain,
do_pml_Lo_MR, amrex::IntVect::TheUnitVector()));
}
}
@@ -195,9 +197,10 @@ WarpX::InitNCICorrector ()
// Initialize Godfrey filters
// Same filter for fields Ex, Ey and Bz
- nci_godfrey_filter_exeybz[lev].reset( new NCIGodfreyFilter(godfrey_coeff_set::Ex_Ey_Bz, cdtodz, WarpX::l_lower_order_in_v) );
+ const bool nodal_gather = (l_lower_order_in_v == 0);
+ nci_godfrey_filter_exeybz[lev].reset( new NCIGodfreyFilter(godfrey_coeff_set::Ex_Ey_Bz, cdtodz, nodal_gather) );
// Same filter for fields Bx, By and Ez
- nci_godfrey_filter_bxbyez[lev].reset( new NCIGodfreyFilter(godfrey_coeff_set::Bx_By_Ez, cdtodz, WarpX::l_lower_order_in_v) );
+ nci_godfrey_filter_bxbyez[lev].reset( new NCIGodfreyFilter(godfrey_coeff_set::Bx_By_Ez, cdtodz, nodal_gather) );
// Compute Godfrey filters stencils
nci_godfrey_filter_exeybz[lev]->ComputeStencils();
nci_godfrey_filter_bxbyez[lev]->ComputeStencils();
@@ -249,9 +252,9 @@ WarpX::InitOpenbc ()
BoxList bl{IndexType::TheNodeType()};
for (int i = 0; i < nprocs; ++i)
{
- bl.push_back(Box(IntVect(alllohi[6*i ],alllohi[6*i+1],alllohi[6*i+2]),
- IntVect(alllohi[6*i+3],alllohi[6*i+4],alllohi[6*i+5]),
- IndexType::TheNodeType()));
+ bl.push_back(Box(IntVect(alllohi[6*i ],alllohi[6*i+1],alllohi[6*i+2]),
+ IntVect(alllohi[6*i+3],alllohi[6*i+4],alllohi[6*i+5]),
+ IndexType::TheNodeType()));
}
BoxArray ba{bl};
@@ -284,13 +287,13 @@ WarpX::InitOpenbc ()
#endif
for (MFIter mfi(phi); mfi.isValid(); ++mfi)
{
- const Box& bx = mfi.validbox();
- warpx_compute_E(bx.loVect(), bx.hiVect(),
- BL_TO_FORTRAN_3D(phi[mfi]),
- BL_TO_FORTRAN_3D((*Efield[lev][0])[mfi]),
- BL_TO_FORTRAN_3D((*Efield[lev][1])[mfi]),
- BL_TO_FORTRAN_3D((*Efield[lev][2])[mfi]),
- dx);
+ const Box& bx = mfi.validbox();
+ warpx_compute_E(bx.loVect(), bx.hiVect(),
+ BL_TO_FORTRAN_3D(phi[mfi]),
+ BL_TO_FORTRAN_3D((*Efield[lev][0])[mfi]),
+ BL_TO_FORTRAN_3D((*Efield[lev][1])[mfi]),
+ BL_TO_FORTRAN_3D((*Efield[lev][2])[mfi]),
+ dx);
}
}
#endif
@@ -299,9 +302,9 @@ void
WarpX::InitLevelData (int lev, Real time)
{
for (int i = 0; i < 3; ++i) {
- current_fp[lev][i]->setVal(0.0);
- Efield_fp[lev][i]->setVal(0.0);
- Bfield_fp[lev][i]->setVal(0.0);
+ current_fp[lev][i]->setVal(0.0);
+ Efield_fp[lev][i]->setVal(0.0);
+ Bfield_fp[lev][i]->setVal(0.0);
}
if (lev > 0) {