pml cubic sigma

RevathiJambunathan · Sep 16, 2023 · e43f9b8 · e43f9b8
1 parent 2b0631d
commit e43f9b8
Show file tree

Hide file tree

Showing 3 changed files with 39 additions and 4 deletions.
diff --git a/Source/BoundaryConditions/PML.cpp b/Source/BoundaryConditions/PML.cpp
@@ -19,6 +19,7 @@
 #include "Utils/WarpXProfilerWrapper.H"
 #include "Utils/Parser/ParserUtils.H"
 #include "WarpX.H"
+#include "Particles/PulsarParameters.H"
 
 #include <ablastr/utils/Communication.H>
 
@@ -71,18 +72,34 @@ namespace
         Real* p_sigma_cumsum = sigma_cumsum.data();
         Real* p_sigma_star = sigma_star.data();
         Real* p_sigma_star_cumsum = sigma_star_cumsum.data();
+#ifdef PULSAR
+        int pml_cubic_sigma = Pulsar::m_pml_cubic_sigma;
+#endif
         amrex::ParallelFor(N, [=] AMREX_GPU_DEVICE (int i) noexcept
         {
             i += olo;
             Real offset = static_cast<Real>(glo-i);
             p_sigma[i-slo] = fac*(offset*offset);
             // sigma_cumsum is the analytical integral of sigma function at same points than sigma
             p_sigma_cumsum[i-slo] = (fac*(offset*offset*offset)/3._rt)/v_sigma;
+#ifdef PULSAR
+            if (pml_cubic_sigma == 1) {
+                p_sigma[i-slo] = fac*(offset*offset*offset);
+                // sigma_cumsum is the analytical integral of sigma function at same points than sigma
+                p_sigma_cumsum[i-slo] = (fac*(offset*offset*offset*offset)/4._rt)/v_sigma;
+            }
+#endif
             if (i <= ohi+1) {
                 offset = static_cast<Real>(glo-i) - 0.5_rt;
                 p_sigma_star[i-sslo] = fac*(offset*offset);
                 // sigma_star_cumsum is the analytical integral of sigma function at same points than sigma_star
                 p_sigma_star_cumsum[i-sslo] = (fac*(offset*offset*offset)/3._rt)/v_sigma;
+#ifdef PULSAR
+                if (pml_cubic_sigma == 1) {
+                    p_sigma_star[i-sslo] = fac*(offset*offset*offset);
+                    p_sigma_star_cumsum[i-sslo] = (fac*(offset*offset*offset*offset)/4._rt)/v_sigma;
+                }
+#endif
             }
         });
     }
@@ -100,16 +117,31 @@ namespace
         Real* p_sigma_cumsum = sigma_cumsum.data();
         Real* p_sigma_star = sigma_star.data();
         Real* p_sigma_star_cumsum = sigma_star_cumsum.data();
+#ifdef PULSAR
+        int pml_cubic_sigma = Pulsar::m_pml_cubic_sigma;
+#endif
         amrex::ParallelFor(N, [=] AMREX_GPU_DEVICE (int i) noexcept
         {
             i += olo;
             Real offset = static_cast<Real>(i-ghi-1);
             p_sigma[i-slo] = fac*(offset*offset);
             p_sigma_cumsum[i-slo] = (fac*(offset*offset*offset)/3._rt)/v_sigma;
+#ifdef PULSAR
+            if (pml_cubic_sigma == 1) {
+                p_sigma[i-slo] = fac*(offset*offset*offset);
+                p_sigma_cumsum[i-slo] = (fac*(offset*offset*offset*offset)/4._rt)/v_sigma;
+            }
+#endif
             if (i <= ohi+1) {
                 offset = static_cast<Real>(i-ghi) - 0.5_rt;
                 p_sigma_star[i-sslo] = fac*(offset*offset);
                 p_sigma_star_cumsum[i-sslo] = (fac*(offset*offset*offset)/3._rt)/v_sigma;
+#ifdef PULSAR
+            if (pml_cubic_sigma == 1) {
+                p_sigma_star[i-sslo] = fac*(offset*offset*offset);
+                p_sigma_star_cumsum[i-sslo] = (fac*(offset*offset*offset*offset)/4._rt)/v_sigma;
+            }
+#endif
             }
         });
     }

diff --git a/Source/Particles/PulsarParameters.H b/Source/Particles/PulsarParameters.H
@@ -701,6 +701,7 @@ public:
     static amrex::Real m_gatherbuffer_max;
     static amrex::Real m_depositbuffer_min;
     static amrex::Real m_depositbuffer_max;
+    static int m_pml_cubic_sigma;
 private:
 };
 

diff --git a/Source/Particles/PulsarParameters.cpp b/Source/Particles/PulsarParameters.cpp
@@ -138,10 +138,11 @@ amrex::Real Pulsar::m_limit_GJfactor = 1.;
 int Pulsar::m_use_BC_smoothening;
 amrex::Real Pulsar::m_min_BC_radius;
 amrex::Real Pulsar::m_BC_width;
-amrex::Real Pulsar::m_gatherbuffer_min;
-amrex::Real Pulsar::m_gatherbuffer_max;
-amrex::Real Pulsar::m_depositbuffer_min;
-amrex::Real Pulsar::m_depositbuffer_max;
+amrex::Real Pulsar::m_gatherbuffer_min = 0.;
+amrex::Real Pulsar::m_gatherbuffer_max = 0.;
+amrex::Real Pulsar::m_depositbuffer_min = 0.;
+amrex::Real Pulsar::m_depositbuffer_max = 0.;
+int Pulsar::m_pml_cubic_sigma;
 
 Pulsar::Pulsar ()
 {
@@ -372,6 +373,7 @@ Pulsar::ReadParameters () {
     pp.query("gatherbuffer_max",m_gatherbuffer_max);
     pp.query("depositbuffer_min",m_depositbuffer_min);
     pp.query("depositbuffer_max",m_depositbuffer_max);
+    pp.get("pml_cubic_sigma", m_pml_cubic_sigma);
 }