Add calculation of beam perveance.

BLAST-ImpactX · Feb 11, 2025 · 1ea8528 · 1ea8528
1 parent 99f69b7
commit 1ea8528
Show file tree

Hide file tree

Showing 2 changed files with 19 additions and 6 deletions.
diff --git a/src/initialization/InitDistribution.cpp b/src/initialization/InitDistribution.cpp
@@ -471,6 +471,14 @@ namespace impactx
         }
         else if (track == "envelope")
         {
+            // For treating 3D space charge in bunched beams (not yet implemented in envelope mode)
+            //amrex::ParticleReal bunch_charge = 0.0;  // Bunch charge (C)
+            //pp_dist.getWithParser("charge", bunch_charge);
+
+            // For treating 2D space charge in unbunched beams
+            amrex::ParticleReal beam_current = 0.0;  // Beam current (A)
+            pp_dist.getWithParser("current", beam_current);
+
             amr_data->track_envelope.m_ref = initialization::read_reference_particle(pp_dist);
             auto dist = initialization::read_distribution(pp_dist);
             amr_data->track_envelope.m_cm = impactx::initialization::create_covariance_matrix(dist);

diff --git a/src/particles/spacecharge/EnvelopeSpaceChargePush.cpp b/src/particles/spacecharge/EnvelopeSpaceChargePush.cpp
@@ -20,7 +20,7 @@ namespace impactx::spacecharge
         void
         envelope_space_charge2D_push (
         [[maybe_unused]] RefPart const & refpart,
-        [[maybe_unused]] Map6x6 & cm,
+        Map6x6 & cm,
         [[maybe_unused]] amrex::ParticleReal current,
         amrex::ParticleReal ds
     )
@@ -30,13 +30,18 @@ namespace impactx::spacecharge
         // initialize the linear transport map
         Map6x6 R = Map6x6::Identity();
 
-        // access reference particle values to find beta*gamma^2
-        //amrex::ParticleReal const pt_ref = refpart.pt;
-        //amrex::ParticleReal const betgam2 = std::pow(pt_ref, 2) - 1.0_prt;
+        // physical constants and reference quantities
+        amrex::ParticleReal const c = ablastr::constant::SI::c;
+        amrex::ParticleReal const ep0 = ablastr::constant::SI::ep0;
+        amrex::ParticleReal const pi = ablastr::constant::math::pi;
+        amrex::ParticleReal const mass = refpart.mass;
+        amrex::ParticleReal const charge = refpart.charge;
+        amrex::ParticleReal const pt_ref = refpart.pt;
+        amrex::ParticleReal const betgam2 = std::pow(pt_ref, 2) - 1.0_prt;
 
         // evaluate the beam space charge perveance from current
-        // TODO
-        amrex::ParticleReal const Kpv=0.0_prt;
+        amrex::ParticleReal const IA = 4.0_prt*pi*ep0*mass*pow(c,3)/charge;
+        amrex::ParticleReal const Kpv = (current/IA) * 2.0_prt/betgam2;
 
         // evaluate the linear transfer map
         amrex::ParticleReal const sigma2 = cm(1,1)*cm(3,3)-cm(1,3)*cm(1,3);