1
0
Fork 0
mirror of https://github.com/gwm17/catima.git synced 2024-11-26 12:08:52 -05:00
This commit is contained in:
hrocho 2018-11-01 01:49:48 +01:00
parent cb8cb30f74
commit 251f1b8c1e
5 changed files with 15 additions and 8 deletions

View File

@ -113,7 +113,7 @@ double bethek_dedx_e(Projectile &p, const Target &t, const Config &c, double I){
double result = (f2)*barkas + LS - delta/2.; double result = (f2)*barkas + LS - delta/2.;
result *=f1; result *=f1;
if( (p.T>50000.0) && (!c.dedx&corrections::no_highenergy)){ if( (p.T>50000.0) && !(c.dedx&corrections::no_highenergy)){
result += pair_production(p,t); result += pair_production(p,t);
result += bremsstrahlung(p,t); result += bremsstrahlung(p,t);
} }
@ -413,8 +413,8 @@ double bethek_lindhard_X(const Projectile &p){
else{ // ultrarelativistic limit else{ // ultrarelativistic limit
} }
return 2*bethek_lindhard(p) - sum - beta2; double res = 2*bethek_lindhard(p) - sum - beta2;
//return sum; return (res>=0)?res:0.0;
} }
double pair_production(const Projectile &p, const Target &t){ double pair_production(const Projectile &p, const Target &t){

View File

@ -297,6 +297,7 @@ class corrections(IntEnum):
no_barkas = 1 no_barkas = 1
no_lindhard = 2 no_lindhard = 2
no_shell_correction = 4 no_shell_correction = 4
no_highenergy = 8
cdef class Config: cdef class Config:
cdef catimac.Config cbase cdef catimac.Config cbase

View File

@ -8,7 +8,7 @@ constexpr double Ezero = 1E-3; // lowest E to calculate, below taken as 0
constexpr double logEmin = -3; // log of minimum energy constexpr double logEmin = -3; // log of minimum energy
constexpr double logEmax = 7.0; // log of max energy constexpr double logEmax = 7.0; // log of max energy
constexpr int max_datapoints = 500; // how many datapoints between logEmin and logEmax constexpr int max_datapoints = 500; // how many datapoints between logEmin and logEmax
constexpr int max_storage_data = 100; // number of datapoints which can be stored in cache constexpr int max_storage_data = 60; // number of datapoints which can be stored in cache
constexpr double numeric_epsilon = std::numeric_limits<double>::epsilon(); constexpr double numeric_epsilon = std::numeric_limits<double>::epsilon();
/// required integration precision (relative units) /// required integration precision (relative units)

View File

@ -53,8 +53,10 @@ namespace catima{
int index(double v)const noexcept{ int index(double v)const noexcept{
double lxval = (log(v/values[0])/M_LN10); double lxval = (log(v/values[0])/M_LN10);
if(v<values[0] || step==0.0)return -1; if(v<values[0] || step==0.0)return -1;
if(v>=values[N-1])return N-1; if(v>=values[N-1]-numeric_epsilon)return N-1;
return static_cast<int> (std::floor(lxval/step)); int i = static_cast<int> (std::floor(lxval/step));
if(v >= values[i+1]-numeric_epsilon)i++; // this is correction for floating point precision
return i;
}; };
std::size_t num; std::size_t num;
}; };
@ -65,7 +67,9 @@ namespace catima{
int EnergyTable_index(const EnergyTable<N> &table, double val){ int EnergyTable_index(const EnergyTable<N> &table, double val){
if(val<table.values[0] || val>table.values[table.num-1])return -1; if(val<table.values[0] || val>table.values[table.num-1])return -1;
double lxval = (log(val/table.values[0])/M_LN10); double lxval = (log(val/table.values[0])/M_LN10);
return static_cast<int>( std::floor(lxval/table.step)); int i = static_cast<int>( std::floor(lxval/table.step));
if(val >= table.values[i+1]-numeric_epsilon)i++; // this is correction for floating point precision
return i;
} }
template<int N> template<int N>

View File

@ -131,11 +131,13 @@ const lest::test specification[] =
EXPECT(EnergyTable_index(catima::energy_table, 0.0)==-1); EXPECT(EnergyTable_index(catima::energy_table, 0.0)==-1);
for(int i:{5,10,100,498}){ for(int i=0;i<catima::max_datapoints-1;i++){
val = catima::energy_table.values[i]; val = catima::energy_table.values[i];
dif = catima::energy_table.values[i+1] - val; dif = catima::energy_table.values[i+1] - val;
EXPECT(EnergyTable_index(catima::energy_table, val)==i); EXPECT(EnergyTable_index(catima::energy_table, val)==i);
EXPECT(EnergyTable_index(catima::energy_table, val+0.5*dif)==i); EXPECT(EnergyTable_index(catima::energy_table, val+0.5*dif)==i);
EXPECT(catima::energy_table.index(val)==i);
EXPECT(catima::energy_table.index(val+0.5*dif)==i);
} }
} }
}; };