Skip to content

Commit

Permalink
max radius for setting uperpB to zero in drifting frame
Browse files Browse the repository at this point in the history
  • Loading branch information
@RevathiJambunathan
RevathiJambunathan committed Nov 6, 2023
1 parent 39352f1 commit 7de20d9
Show file tree
Hide file tree
Showing 7 changed files with 35 additions and 14 deletions.
10 changes: 6 additions & 4 deletions Source/Particles/PhysicalParticleContainer.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2856,6 +2856,7 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
amrex::Real gather_buffer_boxmax = Pulsar::m_gatherbuffer_max;
amrex::Real re_scaledratio = Pulsar::m_re_scaledratio;
int do_zero_uperpB_driftframe = Pulsar::m_do_zero_uperpB_driftframe;
amrex::Real rmax_zero_uperp_driftframe = Pulsar::m_rmax_zero_uperpB_driftframe;
#endif

#ifdef AMREX_USE_OMP
Expand Down Expand Up @@ -3020,7 +3021,7 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
Exp, Eyp, Ezp, Bxp,
Byp, Bzp, qp, m,
#ifdef PULSAR
re_scaledratio, do_zero_uperpB_driftframe,
re_scaledratio, do_zero_uperpB_driftframe, r_p, rmax_zero_uperp_driftframe,
#endif
dt);
} else if (pusher_algo == ParticlePusherAlgo::Boris) {
Expand All @@ -3030,7 +3031,7 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
Exp, Eyp, Ezp, Bxp,
Byp, Bzp, qp, m,
#ifdef PULSAR
do_zero_uperpB_driftframe,
do_zero_uperpB_driftframe, r_p, rmax_zero_uperp_driftframe,
#endif
dt);
} else if (pusher_algo == ParticlePusherAlgo::Vay) {
Expand Down Expand Up @@ -3170,6 +3171,7 @@ PhysicalParticleContainer::PushPX (WarpXParIter& pti,
amrex::Real gather_buffer_boxmax = Pulsar::m_gatherbuffer_max;
amrex::Real re_scaledratio = Pulsar::m_re_scaledratio;
int do_zero_uperpB_driftframe = Pulsar::m_do_zero_uperpB_driftframe;
amrex::Real rmax_zero_uperp_driftframe = Pulsar::m_rmax_zero_uperpB_driftframe;
#endif


Expand Down Expand Up @@ -3421,7 +3423,7 @@ PhysicalParticleContainer::PushPX (WarpXParIter& pti,
t_chi_max,
#endif
#ifdef PULSAR
re_scaledratio, do_zero_uperpB_driftframe,
re_scaledratio, do_zero_uperpB_driftframe, r_p, rmax_zero_uperp_driftframe,
#endif
dt);
// PulsarPartDiagData);
Expand Down Expand Up @@ -3461,7 +3463,7 @@ PhysicalParticleContainer::PushPX (WarpXParIter& pti,
m, q, pusher_algo, do_crr, do_copy,
t_chi_max,
#ifdef PULSAR
re_scaledratio, do_zero_uperpB_driftframe,
re_scaledratio, do_zero_uperpB_driftframe, r_p, rmax_zero_uperp_driftframe,
#endif
dt);
}
Expand Down
1 change: 1 addition & 0 deletions Source/Particles/PulsarParameters.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -714,6 +714,7 @@ public:
static amrex::Real m_beta_rec_RR;
static amrex::Real m_re_scaledratio;
static int m_do_zero_uperpB_driftframe;
static amrex::Real m_rmax_zero_uperpB_driftframe;
private:
};

Expand Down
2 changes: 2 additions & 0 deletions Source/Particles/PulsarParameters.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -154,6 +154,7 @@ amrex::Real Pulsar::m_BLC;
amrex::Real Pulsar::m_beta_rec_RR;
amrex::Real Pulsar::m_re_scaledratio = 1.;
int Pulsar::m_do_zero_uperpB_driftframe = 0;
amrex::Real Pulsar::m_rmax_zero_uperpB_driftframe;

Pulsar::Pulsar ()
{
Expand Down Expand Up @@ -406,6 +407,7 @@ Pulsar::ReadParameters () {
}
pp.get("do_zero_uperpB_driftframe",m_do_zero_uperpB_driftframe);
amrex::Print() << " do zero uperp B " << m_do_zero_uperpB_driftframe << "\n";
pp.get("rmax_zero_uperpB_driftframe",m_rmax_zero_uperpB_driftframe);
}


Expand Down
14 changes: 10 additions & 4 deletions Source/Particles/Pusher/PushSelector.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -68,6 +68,8 @@ void doParticlePush(const GetParticlePosition& GetPosition,
#ifdef PULSAR
const amrex::Real re_scaledratio,
const int do_zero_uperpB_driftframe,
const amrex::ParticleReal r_p,
const amrex::Real rmax_zero_uperp_driftframe,
#endif
const amrex::Real dt, amrex::Real * PulsarDiag = nullptr)
{
Expand All @@ -87,7 +89,8 @@ void doParticlePush(const GetParticlePosition& GetPosition,
Ex, Ey, Ez, Bx,
By, Bz, qp, m,
#ifdef PULSAR
re_scaledratio, do_zero_uperpB_driftframe,
re_scaledratio, do_zero_uperpB_driftframe, r_p,
rmax_zero_uperp_driftframe,
#endif
dt);
}
Expand All @@ -96,7 +99,8 @@ void doParticlePush(const GetParticlePosition& GetPosition,
Ex, Ey, Ez, Bx,
By, Bz, qp, m,
#ifdef PULSAR
do_zero_uperpB_driftframe,
do_zero_uperpB_driftframe, r_p,
rmax_zero_uperp_driftframe,
#endif
dt);
}
Expand All @@ -112,7 +116,8 @@ void doParticlePush(const GetParticlePosition& GetPosition,
Ex, Ey, Ez, Bx,
By, Bz, qp, m,
#ifdef PULSAR
re_scaledratio, do_zero_uperpB_driftframe,
re_scaledratio, do_zero_uperpB_driftframe, r_p,
rmax_zero_uperp_driftframe,
#endif
dt);
amrex::ParticleReal x, y, z;
Expand All @@ -125,7 +130,8 @@ void doParticlePush(const GetParticlePosition& GetPosition,
Ex, Ey, Ez, Bx,
By, Bz, qp, m,
#ifdef PULSAR
do_zero_uperpB_driftframe,
do_zero_uperpB_driftframe, r_p,
rmax_zero_uperp_driftframe,
#endif
dt);
amrex::ParticleReal x, y, z;
Expand Down
8 changes: 5 additions & 3 deletions Source/Particles/Pusher/UpdateMomentumBoris.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,8 @@ void UpdateMomentumBoris(
const amrex::ParticleReal q, const amrex::ParticleReal m,
#ifdef PULSAR
const int do_zero_uperpB_driftframe,
const amrex::ParticleReal r_p,
const amrex::Real rmax_zero_uperpB_driftframe,
#endif
const amrex::Real dt, amrex::Real * PulsarDiag = nullptr )
{
Expand Down Expand Up @@ -68,7 +70,7 @@ void UpdateMomentumBoris(


#ifdef PULSAR
if (do_zero_uperpB_driftframe) {
if (do_zero_uperpB_driftframe && (r_p < rmax_zero_uperpB_driftframe) ) {
// 1(a) Find drift velocity
amrex::Real B_sq = Bx*Bx + By*By + Bz*Bz;

Expand All @@ -83,14 +85,14 @@ void UpdateMomentumBoris(
+ std::sqrt( ( (B_sq*c2) - E_sq) * ( (B_sq*c2) - E_sq ) + 4._prt * EdotB_sq * c2 ) );

// only if E << B -- else fix with thesis
if ( (B_sq + Eprime_sq/c2) > 0 ) {
if ( (B_sq + Eprime_sq/c2) > 1.e-30 ) {
amrex::Real vx_drift = (Ey*Bz - Ez*By)/(B_sq + Eprime_sq/c2);
amrex::Real vy_drift = (Ez*Bx - Ex*Bz)/(B_sq + Eprime_sq/c2);
amrex::Real vz_drift = (Ex*By - Ey*Bx)/(B_sq + Eprime_sq/c2);

// 1(b) unit vector along drift velocity
amrex::Real drift_mag = std::sqrt(vx_drift*vx_drift + vy_drift*vy_drift + vz_drift*vz_drift);
if (drift_mag < PhysConst::c) {
if (drift_mag < (0.99 * PhysConst::c) ) {
amrex::Real vxd_hat = vx_drift/drift_mag;
amrex::Real vyd_hat = vy_drift/drift_mag;
amrex::Real vzd_hat = vz_drift/drift_mag;
Expand Down
4 changes: 3 additions & 1 deletion Source/Particles/Pusher/UpdateMomentumBorisWithRadiationReaction.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,8 @@ void UpdateMomentumBorisWithRadiationReaction(
#ifdef PULSAR
const amrex::Real re_scaledratio,
const int do_zero_uperpB_driftframe,
const amrex::ParticleReal r_p,
const amrex::Real rmax_zero_uperpB_driftframe,
#endif
const amrex::Real dt )
{
Expand All @@ -49,7 +51,7 @@ void UpdateMomentumBorisWithRadiationReaction(
Bx, By, Bz,
q, m,
#ifdef PULSAR
do_zero_uperpB_driftframe,
do_zero_uperpB_driftframe, r_p, rmax_zero_uperpB_driftframe,
#endif
dt );

Expand Down
10 changes: 8 additions & 2 deletions Source/Particles/RigidInjectedParticleContainer.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -411,6 +411,7 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt,
#ifdef PULSAR
amrex::Real re_scaledratio = Pulsar::m_re_scaledratio;
int do_zero_uperpB_driftframe = Pulsar::m_do_zero_uperpB_driftframe;
amrex::Real rmax_zero_uperp_driftframe = Pulsar::m_rmax_zero_uperpB_driftframe;
#endif

amrex::ParallelFor( np, [=] AMREX_GPU_DEVICE (long ip)
Expand All @@ -422,6 +423,9 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt,
amrex::ParticleReal xp, yp, zp;
getPosition(ip, xp, yp, zp);

#ifdef PULSAR
amrex::ParticleReal r_p = std::sqrt(xp*xp + yp*yp + zp*zp);
#endif
amrex::ParticleReal Exp = 0._prt, Eyp = 0._prt, Ezp = 0._prt;
amrex::ParticleReal Bxp = 0._prt, Byp = 0._prt, Bzp = 0._prt;

Expand All @@ -442,15 +446,17 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt,
Byp, Bzp, qp, m,
#ifdef PULSAR
re_scaledratio,
do_zero_uperpB_driftframe,
do_zero_uperpB_driftframe, r_p,
rmax_zero_uperp_driftframe,
#endif
dt);
} else if (pusher_algo == ParticlePusherAlgo::Boris) {
UpdateMomentumBoris( uxpp[ip], uypp[ip], uzpp[ip],
Exp, Eyp, Ezp, Bxp,
Byp, Bzp, qp, m,
#ifdef PULSAR
do_zero_uperpB_driftframe,
do_zero_uperpB_driftframe, r_p,
rmax_zero_uperp_driftframe,
#endif
dt);
} else if (pusher_algo == ParticlePusherAlgo::Vay) {
Expand Down

0 comments on commit 7de20d9

Please sign in to comment.