Merge pull request #461 from ChrisArran/ChrisArran-photonEnergyExtrap…

…olation Bug fix: photon energy extrapolation
Warwick-Plasma · Mar 13, 2023 · 05fdd6f · 05fdd6f
2 parents f371258 + ac93b12
commit 05fdd6f
Show file tree

Hide file tree

Showing 3 changed files with 30 additions and 9 deletions.
diff --git a/epoch1d/src/physics_packages/photons.F90 b/epoch1d/src/physics_packages/photons.F90
@@ -916,10 +916,17 @@ FUNCTION calculate_photon_energy(rand_seed, eta, generating_gamma)
 
     REAL(num) :: calculate_photon_energy
     REAL(num), INTENT(IN) :: rand_seed, eta, generating_gamma
-    REAL(num) :: chi_final
+    REAL(num) :: eta_min, chi_tmp, chi_final
 
-    chi_final = find_value_from_table_alt(eta, rand_seed, &
-        n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+    eta_min = 10.0_num**MINVAL(log_eta)
+    IF (eta < eta_min) THEN ! Extrapolate downwards with chi \propto eta^2
+      chi_tmp = find_value_from_table_alt(eta_min, rand_seed, &
+          n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+      chi_final = chi_tmp * (eta / eta_min)**2
+    ELSE
+      chi_final = find_value_from_table_alt(eta, rand_seed, &
+          n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+    ENDIF
 
     calculate_photon_energy = (2.0_num * chi_final / eta) * generating_gamma &
         * m0 * c**2

diff --git a/epoch2d/src/physics_packages/photons.F90 b/epoch2d/src/physics_packages/photons.F90
@@ -929,10 +929,17 @@ FUNCTION calculate_photon_energy(rand_seed, eta, generating_gamma)
 
     REAL(num) :: calculate_photon_energy
     REAL(num), INTENT(IN) :: rand_seed, eta, generating_gamma
-    REAL(num) :: chi_final
+    REAL(num) :: eta_min, chi_tmp, chi_final  
 
-    chi_final = find_value_from_table_alt(eta, rand_seed, &
-        n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+    eta_min = 10.0_num**MINVAL(log_eta)
+    IF (eta < eta_min) THEN ! Extrapolate downwards with chi \propto eta^2
+      chi_tmp = find_value_from_table_alt(eta_min, rand_seed, &
+          n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+      chi_final = chi_tmp * (eta / eta_min)**2
+    ELSE
+      chi_final = find_value_from_table_alt(eta, rand_seed, &
+          n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+    ENDIF
 
     calculate_photon_energy = (2.0_num * chi_final / eta) * generating_gamma &
         * m0 * c**2

diff --git a/epoch3d/src/physics_packages/photons.F90 b/epoch3d/src/physics_packages/photons.F90
@@ -942,10 +942,17 @@ FUNCTION calculate_photon_energy(rand_seed, eta, generating_gamma)
 
     REAL(num) :: calculate_photon_energy
     REAL(num), INTENT(IN) :: rand_seed, eta, generating_gamma
-    REAL(num) :: chi_final
+    REAL(num) :: eta_min, chi_tmp, chi_final  
 
-    chi_final = find_value_from_table_alt(eta, rand_seed, &
-        n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+    eta_min = 10.0_num**MINVAL(log_eta)
+    IF (eta < eta_min) THEN ! Extrapolate downwards with chi \propto eta^2
+      chi_tmp = find_value_from_table_alt(eta_min, rand_seed, &
+          n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+      chi_final = chi_tmp * (eta / eta_min)**2
+    ELSE
+      chi_final = find_value_from_table_alt(eta, rand_seed, &
+          n_sample_eta, n_sample_chi, log_eta, log_chi, p_photon_energy)
+    ENDIF
 
     calculate_photon_energy = (2.0_num * chi_final / eta) * generating_gamma &
         * m0 * c**2