Skip to content

Commit

Permalink
scale re for radiation using gamma_rad
Browse files Browse the repository at this point in the history
  • Loading branch information
@RevathiJambunathan
RevathiJambunathan committed Oct 20, 2023
1 parent 813ba90 commit 1f0233f
Show file tree
Hide file tree
Showing 6 changed files with 42 additions and 21 deletions.
12 changes: 5 additions & 7 deletions Source/Particles/PhysicalParticleContainer.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2854,8 +2854,7 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
amrex::Real BC_width = Pulsar::m_BC_width;
amrex::Real gather_buffer_boxmin = Pulsar::m_gatherbuffer_min;
amrex::Real gather_buffer_boxmax = Pulsar::m_gatherbuffer_max;
amrex::Real gammarad_real = Pulsar::m_gammarad_real;
amrex::Real gammarad_scaled = Pulsar::m_gammarad_scaled;
amrex::Real re_scaledratio = Pulsar::m_re_scaledratio;
#endif

#ifdef AMREX_USE_OMP
Expand Down Expand Up @@ -3020,7 +3019,7 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
Exp, Eyp, Ezp, Bxp,
Byp, Bzp, qp, m,
#ifdef PULSAR
gammarad_real, gammarad_scaled,
re_scaledratio,
#endif
dt);
} else if (pusher_algo == ParticlePusherAlgo::Boris) {
Expand Down Expand Up @@ -3164,8 +3163,7 @@ PhysicalParticleContainer::PushPX (WarpXParIter& pti,
amrex::Real BC_width = Pulsar::m_BC_width;
amrex::Real gather_buffer_boxmin = Pulsar::m_gatherbuffer_min;
amrex::Real gather_buffer_boxmax = Pulsar::m_gatherbuffer_max;
amrex::Real gammarad_real = Pulsar::m_gammarad_real;
amrex::Real gammarad_scaled = Pulsar::m_gammarad_scaled;
amrex::Real re_scaledratio = Pulsar::m_re_scaledratio;
#endif


Expand Down Expand Up @@ -3417,7 +3415,7 @@ PhysicalParticleContainer::PushPX (WarpXParIter& pti,
t_chi_max,
#endif
#ifdef PULSAR
gammarad_real, gammarad_scaled,
re_scaledratio,
#endif
dt);
// PulsarPartDiagData);
Expand Down Expand Up @@ -3457,7 +3455,7 @@ PhysicalParticleContainer::PushPX (WarpXParIter& pti,
m, q, pusher_algo, do_crr, do_copy,
t_chi_max,
#ifdef PULSAR
gammarad_real, gammarad_scaled,
re_scaledratio,
#endif
dt);
}
Expand Down
6 changes: 6 additions & 0 deletions Source/Particles/PulsarParameters.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -707,6 +707,12 @@ public:
static amrex::Real m_gammarad_scaled;
static amrex::Real m_damping_strength;
static amrex::Real m_totalcells_injectionring;
static int m_do_scale_re_RR;
static amrex::Real m_re_scaled;
static amrex::Real m_gammarad_RR;
static amrex::Real m_BLC;
static amrex::Real m_beta_rec_RR;
static amrex::Real m_re_scaledratio;
private:
};

Expand Down
19 changes: 19 additions & 0 deletions Source/Particles/PulsarParameters.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -147,6 +147,12 @@ amrex::Real Pulsar::m_gammarad_real = 1.;
amrex::Real Pulsar::m_gammarad_scaled = 1.;
amrex::Real Pulsar::m_damping_strength = 4.;
amrex::Real Pulsar::m_totalcells_injectionring = 0;
int Pulsar::m_do_scale_re_RR = 0;
amrex::Real Pulsar::m_re_scaled;
amrex::Real Pulsar::m_gammarad_RR;
amrex::Real Pulsar::m_BLC;
amrex::Real Pulsar::m_beta_rec_RR;
amrex::Real Pulsar::m_re_scaledratio = 1.;

Pulsar::Pulsar ()
{
Expand Down Expand Up @@ -383,6 +389,19 @@ Pulsar::ReadParameters () {
pp.get("damping_strength", m_damping_strength);
amrex::Print() << " damping strength in PML " << m_damping_strength << "\n";
amrex::Print() << " cubic pml flag : " << m_pml_cubic_sigma << "\n";
// for radiation reaction with scaled electron radius
pp.get("scale_re_for_RR", m_do_scale_re_RR);
if (m_do_scale_re_RR == 1) {
pp.get("gammarad_RR",m_gammarad_RR); // gammathreshold at which re is scaled
pp.get("beta_rec_RR",m_beta_rec_RR); // reconnection rate
amrex::Real S = m_RLC/m_R_star; // Scale separation ratio
m_BLC = m_B_star/(S*S*S); // Bfield at LC
m_re_scaled = m_beta_rec_RR * (3./2.) * PhysConst::m_e * PhysConst::c
/ (m_gammarad_RR * m_gammarad_RR * PhysConst::q_e * m_BLC);
amrex::Print() << " re_scaled at gammarad : " << m_gammarad_RR << " is : " << m_re_scaled << "\n";
m_re_scaledratio = m_re_scaled / PhysConst::r_e;
amrex::Print() << " ratio : re_scaled/re : " << m_re_scaledratio << "\n";
}
}


Expand Down
6 changes: 3 additions & 3 deletions Source/Particles/Pusher/PushSelector.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -66,7 +66,7 @@ void doParticlePush(const GetParticlePosition& GetPosition,
const amrex::Real t_chi_max,
#endif
#ifdef PULSAR
const amrex::Real gammarad_real, const amrex::Real gammarad_scaled,
const amrex::Real re_scaledratio,
#endif
const amrex::Real dt, amrex::Real * PulsarDiag = nullptr)
{
Expand All @@ -86,7 +86,7 @@ void doParticlePush(const GetParticlePosition& GetPosition,
Ex, Ey, Ez, Bx,
By, Bz, qp, m,
#ifdef PULSAR
gammarad_real, gammarad_scaled,
re_scaledratio,
#endif
dt);
}
Expand All @@ -107,7 +107,7 @@ void doParticlePush(const GetParticlePosition& GetPosition,
Ex, Ey, Ez, Bx,
By, Bz, qp, m,
#ifdef PULSAR
gammarad_real, gammarad_scaled,
re_scaledratio,
#endif
dt);
amrex::ParticleReal x, y, z;
Expand Down
15 changes: 7 additions & 8 deletions Source/Particles/Pusher/UpdateMomentumBorisWithRadiationReaction.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@ void UpdateMomentumBorisWithRadiationReaction(
const amrex::ParticleReal Bx, const amrex::ParticleReal By, const amrex::ParticleReal Bz,
const amrex::ParticleReal q, const amrex::ParticleReal m,
#ifdef PULSAR
const amrex::Real gammarad_real, const amrex::Real gammarad_scaled,
const amrex::Real re_scaledratio,
#endif
const amrex::Real dt )
{
Expand Down Expand Up @@ -75,18 +75,17 @@ void UpdateMomentumBorisWithRadiationReaction(
//Calculation of auxiliary quantities
const amrex::ParticleReal bdotE = (bx_n*Ex + by_n*Ey + bz_n*Ez);
const amrex::ParticleReal bdotE2 = bdotE*bdotE;
#ifdef PULSAR
// To account for radiation reaction in the scaled pulsar simulation
// we scale the ultra-relativistic term, gamma_n with gamma_real/gamma_scaled
const amrex::ParticleReal gamma_scale_fac = gammarad_scaled/gammarad_real;
const amrex::ParticleReal coeff = gamma_scale_fac*gamma_scale_fac*gamma_n*gamma_n*(fl_q2-bdotE2);
#else
const amrex::ParticleReal coeff = gamma_n*gamma_n*(fl_q2-bdotE2);
#endif

//Radiation reaction constant
const amrex::ParticleReal q_over_mc = q/(m*PhysConst::c);
#ifdef PULSAR
// For a scaled down system, we scale up radiation reaction drag force by upscaling re
// re_scaled is derived by equating the drag force and accelerating force, for a set gammarad < sigma_LC
const amrex::ParticleReal RRcoeff = (2.0_prt/3.0_prt)*PhysConst::r_e*re_scaledratio*q_over_mc*q_over_mc;
#else
const amrex::ParticleReal RRcoeff = (2.0_prt/3.0_prt)*PhysConst::r_e*q_over_mc*q_over_mc;
#endif

//Compute the components of the RR force
const amrex::ParticleReal frx =
Expand Down
5 changes: 2 additions & 3 deletions Source/Particles/RigidInjectedParticleContainer.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -409,8 +409,7 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt,
const auto do_crr = do_classical_radiation_reaction;

#ifdef PULSAR
amrex::Real gammarad_real = Pulsar::m_gammarad_real;
amrex::Real gammarad_scaled = Pulsar::m_gammarad_scaled;
amrex::Real re_scaledratio = Pulsar::m_re_scaledratio;
#endif

amrex::ParallelFor( np, [=] AMREX_GPU_DEVICE (long ip)
Expand Down Expand Up @@ -441,7 +440,7 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt,
Exp, Eyp, Ezp, Bxp,
Byp, Bzp, qp, m,
#ifdef PULSAR
gammarad_real, gammarad_scaled,
re_scaledratio,
#endif
dt);
} else if (pusher_algo == ParticlePusherAlgo::Boris) {
Expand Down

0 comments on commit 1f0233f

Please sign in to comment.