1
0
Fork 0
mirror of https://github.com/gwm17/catima.git synced 2024-11-22 10:18:50 -05:00
This commit is contained in:
hrocho 2019-10-08 20:49:52 +02:00
parent e80e0e26f6
commit e2f67686ea
8 changed files with 470 additions and 1873 deletions

File diff suppressed because it is too large Load Diff

View File

@ -1,73 +1,64 @@
#include "lest.hpp"
#define DOCTEST_CONFIG_IMPLEMENT_WITH_MAIN
#define DOCTEST_CONFIG_SUPER_FAST_ASSERTS
#include "doctest.h"
#include <math.h>
#include "catima/abundance_database.h"
#include "testutils.h"
using namespace std;
using catima::approx;
using namespace abundance;
const lest::test specification[] =
{
CASE("isotope abundaces"){
EXPECT(get_isotopes_num(1) == 2);
EXPECT(get_isotopes_num(27) == 1);
EXPECT(get_isotopes_num(0) == 0);
EXPECT(get_isotopes_num(111) == 0);
TEST_CASE("isotope abundaces"){
CHECK(get_isotopes_num(1) == 2);
CHECK(get_isotopes_num(27) == 1);
CHECK(get_isotopes_num(0) == 0);
CHECK(get_isotopes_num(111) == 0);
EXPECT(get_isotope_a(0,0) == 0);
EXPECT(get_isotope_a(0,1) == 0);
EXPECT(get_isotope_a(1,0) == 1);
EXPECT(get_isotope_a(1,1) == 2);
EXPECT(get_isotope_a(1,2) == 0);
EXPECT(get_isotope_a(14,0) == 28);
EXPECT(get_isotope_a(14,1) == 29);
EXPECT(get_isotope_a(14,2) == 30);
EXPECT(get_isotope_a(14,3) == 0);
EXPECT(get_isotope_a(68,0) == 166);
EXPECT(get_isotope_a(68,1) == 168);
EXPECT(get_isotope_a(120,0) == 0);
CHECK(get_isotope_a(0,0) == 0);
CHECK(get_isotope_a(0,1) == 0);
CHECK(get_isotope_a(1,0) == 1);
CHECK(get_isotope_a(1,1) == 2);
CHECK(get_isotope_a(1,2) == 0);
CHECK(get_isotope_a(14,0) == 28);
CHECK(get_isotope_a(14,1) == 29);
CHECK(get_isotope_a(14,2) == 30);
CHECK(get_isotope_a(14,3) == 0);
CHECK(get_isotope_a(68,0) == 166);
CHECK(get_isotope_a(68,1) == 168);
CHECK(get_isotope_a(120,0) == 0);
EXPECT(get_abundance(1,0) == approx(0.999,0.01));
EXPECT(get_abundance(1,2) == 0.0);
EXPECT(get_abundance(14,0) == approx(0.922,0.01));
EXPECT(get_abundance(14,1) == approx(0.046,0.01));
EXPECT(get_abundance(14,2) == approx(0.0308,0.01));
EXPECT(get_abundance(14,3) == 0.0);
EXPECT(get_abundance(68,0) == approx(0.336,0.01));
EXPECT(get_abundance(68,1) == approx(0.267,0.01));
EXPECT(get_abundance(120,0) == 0);
CHECK(get_abundance(1,0) == approx(0.999,0.01));
CHECK(get_abundance(1,2) == 0.0);
CHECK(get_abundance(14,0) == approx(0.922,0.01));
CHECK(get_abundance(14,1) == approx(0.046,0.01));
CHECK(get_abundance(14,2) == approx(0.0308,0.01));
CHECK(get_abundance(14,3) == 0.0);
CHECK(get_abundance(68,0) == approx(0.336,0.01));
CHECK(get_abundance(68,1) == approx(0.267,0.01));
CHECK(get_abundance(120,0) == 0);
EXPECT(get_isotope(0,0).first == 0);
EXPECT(get_isotope(0,1).first == 0);
EXPECT(get_isotope(1,0).first == 1);
EXPECT(get_isotope(1,1).first == 2);
EXPECT(get_isotope(1,2).first == 0);
EXPECT(get_isotope(14,0).first == 28);
EXPECT(get_isotope(14,1).first == 29);
EXPECT(get_isotope(14,2).first == 30);
EXPECT(get_isotope(14,3).first == 0);
EXPECT(get_isotope(68,0).first == 166);
EXPECT(get_isotope(68,1).first == 168);
EXPECT(get_isotope(120,0).first == 0);
CHECK(get_isotope(0,0).first == 0);
CHECK(get_isotope(0,1).first == 0);
CHECK(get_isotope(1,0).first == 1);
CHECK(get_isotope(1,1).first == 2);
CHECK(get_isotope(1,2).first == 0);
CHECK(get_isotope(14,0).first == 28);
CHECK(get_isotope(14,1).first == 29);
CHECK(get_isotope(14,2).first == 30);
CHECK(get_isotope(14,3).first == 0);
CHECK(get_isotope(68,0).first == 166);
CHECK(get_isotope(68,1).first == 168);
CHECK(get_isotope(120,0).first == 0);
EXPECT(get_isotope(1,0).second == approx(0.999,0.01));
EXPECT(get_isotope(1,2).second == 0.0);
EXPECT(get_isotope(14,0).second == approx(0.922,0.01));
EXPECT(get_isotope(14,1).second == approx(0.046,0.01));
EXPECT(get_isotope(14,2).second == approx(0.0308,0.01));
EXPECT(get_isotope(14,3).second == 0.0);
EXPECT(get_isotope(68,0).second == approx(0.336,0.01));
EXPECT(get_isotope(68,1).second == approx(0.267,0.01));
EXPECT(get_isotope(120,0).second == 0);
CHECK(get_isotope(1,0).second == approx(0.999,0.01));
CHECK(get_isotope(1,2).second == 0.0);
CHECK(get_isotope(14,0).second == approx(0.922,0.01));
CHECK(get_isotope(14,1).second == approx(0.046,0.01));
CHECK(get_isotope(14,2).second == approx(0.0308,0.01));
CHECK(get_isotope(14,3).second == 0.0);
CHECK(get_isotope(68,0).second == approx(0.336,0.01));
CHECK(get_isotope(68,1).second == approx(0.267,0.01));
CHECK(get_isotope(120,0).second == 0);
}
};
int main( int argc, char * argv[] )
{
return lest::run( specification, argc, argv );
}

View File

@ -1,161 +1,156 @@
#include "lest.hpp"
#define DOCTEST_CONFIG_IMPLEMENT_WITH_MAIN
#define DOCTEST_CONFIG_SUPER_FAST_ASSERTS
#include "doctest.h"
#include "testutils.h"
#include <math.h>
#include "catima/catima.h"
#include "catima/calculations.h"
using namespace std;
using catima::approx;
const lest::test specification[] =
{
CASE("nuclear stopping power"){
TEST_CASE("nuclear stopping power"){
catima::Target carbon{12.0107,6};
catima::Projectile p{4.00151,2,2,1};
double dif;
p.T = 0.1/p.A; //0.1MeV
dif = catima::dedx_n(p,carbon) - 14.27;
EXPECT( fabs(dif)< 1);
CHECK( fabs(dif)< 1);
p.T = 1/p.A; //1MeV
dif = catima::dedx_n(p,carbon) - 2.161;
EXPECT( fabs(dif)< 0.1);
CHECK( fabs(dif)< 0.1);
p.T = 10/p.A; //10MeV
dif = catima::dedx_n(p,carbon) - 0.2874;
EXPECT( fabs(dif) < 0.01);
CHECK( fabs(dif) < 0.01);
p.T = 100/p.A; //100MeV
dif = catima::dedx_n(p,carbon) - 0.03455;
EXPECT( fabs(dif) < 0.001);
},
CHECK( fabs(dif) < 0.001);
}
CASE("proton stopping power from srim"){
TEST_CASE("proton stopping power from srim"){
catima::Projectile p{1,1,1,1};
auto he = catima::get_material(2);
auto carbon = catima::get_material(6);
p.T = 1;
EXPECT( catima::sezi_dedx_e(p,he) == approx(283,1));
CHECK( catima::sezi_dedx_e(p,he) == approx(283,1));
p.T = 10;
EXPECT( catima::sezi_dedx_e(p,he) == approx(45.6,1));
CHECK( catima::sezi_dedx_e(p,he) == approx(45.6,1));
p.T = 30;
EXPECT( catima::sezi_dedx_e(p,he) == approx(18.38,1));
CHECK( catima::sezi_dedx_e(p,he) == approx(18.38,1));
p.T = 1;
EXPECT( catima::sezi_dedx_e(p,carbon) == approx(229.5,1));
CHECK( catima::sezi_dedx_e(p,carbon) == approx(229.5,1));
p.T = 10;
EXPECT( catima::sezi_dedx_e(p,carbon) == approx(40.8,1));
CHECK( catima::sezi_dedx_e(p,carbon) == approx(40.8,1));
p.T = 30;
EXPECT( catima::sezi_dedx_e(p,carbon) == approx(16.8,1));
},
CASE("dedx, low energy, from sezi"){
CHECK( catima::sezi_dedx_e(p,carbon) == approx(16.8,1));
}
TEST_CASE("dedx, low energy, from sezi"){
catima::Projectile p{4,2,2,1};
auto carbon = catima::get_material(6);
// He projectile case
// He projectile TEST_CASE
p.T = 1;
EXPECT( catima::sezi_dedx_e(p,carbon)+catima::dedx_n(p,carbon) == approx(922.06).R(0.0001) );
CHECK( catima::sezi_dedx_e(p,carbon)+catima::dedx_n(p,carbon) == approx(922.06).R(0.0001) );
p.T = 3;
EXPECT( catima::sezi_dedx_e(p,carbon)+catima::dedx_n(p,carbon) == approx(433.09).R(0.0001) );
CHECK( catima::sezi_dedx_e(p,carbon)+catima::dedx_n(p,carbon) == approx(433.09).R(0.0001) );
// C projectile case
// C projectile TEST_CASE
p.A = 12;
p.Z = 6;
p.T = 1;
EXPECT( catima::sezi_dedx_e(p,carbon)+catima::dedx_n(p,carbon) == approx( 5792.52).R(0.0001) );
CHECK( catima::sezi_dedx_e(p,carbon)+catima::dedx_n(p,carbon) == approx( 5792.52).R(0.0001) );
p.T = 9.9;
EXPECT( catima::sezi_dedx_e(p,carbon)+catima::dedx_n(p,carbon) == approx(1485.36).R(0.0001) );
CHECK( catima::sezi_dedx_e(p,carbon)+catima::dedx_n(p,carbon) == approx(1485.36).R(0.0001) );
},
CASE("LS check: deltaL values"){
}
TEST_CASE("LS check: deltaL values"){
catima::Projectile p{238,92,92,1};
p.T = 93.1494;
EXPECT(catima::bethek_lindhard(p)== approx(-0.5688,0.0001));
CHECK(catima::bethek_lindhard(p)== approx(-0.5688,0.0001));
p.T = 380.9932;
EXPECT(catima::bethek_lindhard(p)== approx(0.549199,0.0001));
CHECK(catima::bethek_lindhard(p)== approx(0.549199,0.0001));
p.T = 995.368987;
EXPECT(catima::bethek_lindhard(p)== approx(1.106649).R(0.001) );
CHECK(catima::bethek_lindhard(p)== approx(1.106649).R(0.001) );
p.T = 2640.032566;
EXPECT(catima::bethek_lindhard(p)== approx(1.35314).R(0.001) );
CHECK(catima::bethek_lindhard(p)== approx(1.35314).R(0.001) );
p.T = 6091.392448;
EXPECT(catima::bethek_lindhard(p)== approx(1.365643).R(0.001) );
CHECK(catima::bethek_lindhard(p)== approx(1.365643).R(0.001) );
p.T = 37277.695445;
EXPECT(catima::bethek_lindhard(p)== approx(0.689662).R(0.001) );
},
CHECK(catima::bethek_lindhard(p)== approx(0.689662).R(0.001) );
}
CASE("LS check: straggling values"){
TEST_CASE("LS check: straggling values"){
catima::Projectile p{238,92,92,1};
auto f = [&](){return catima::bethek_lindhard_X(p);};
p.T = 93.1494;
EXPECT( f() == approx(1.56898).R(0.01) );
CHECK( f() == approx(1.56898).R(0.01) );
p.T = 380.9932;
EXPECT( f() == approx(1.836008).R(0.01) );
CHECK( f() == approx(1.836008).R(0.01) );
p.T = 996.9855;
EXPECT( f() == approx(1.836528).R(0.01) );
CHECK( f() == approx(1.836528).R(0.01) );
p.T = 2794.4822;
EXPECT( f()== approx(1.768018).R(0.01) );
CHECK( f()== approx(1.768018).R(0.01) );
for(double e:{2000,20000,200000, 9000000, 50000000})
EXPECT(catima::precalculated_lindhard_X(p(e)) >= 0.0);
},
CASE("ultrarelativistic corrections"){
CHECK(catima::precalculated_lindhard_X(p(e)) >= 0.0);
}
TEST_CASE("ultrarelativistic corrections"){
catima::Projectile p{238,92};
catima::Target t{27,13};
EXPECT(catima::pair_production(p(1e3),t) == approx(0.0,1e-3));
EXPECT(catima::bremsstrahlung(p(1e3),t) == approx(0.0,1e-3));
CHECK(catima::pair_production(p(1e3),t) == approx(0.0,1e-3));
CHECK(catima::bremsstrahlung(p(1e3),t) == approx(0.0,1e-3));
EXPECT(catima::pair_production(p(1e6),t) == approx(1900,300));
EXPECT(catima::bremsstrahlung(p(1e6),t) == approx(170,20));
EXPECT(catima::pair_production(p(7e6),t) == approx(19000,5000));
EXPECT(catima::bremsstrahlung(p(7e6),t) == approx(6000,500));
},
CASE("dEdx for compounds"){
CHECK(catima::pair_production(p(1e6),t) == approx(1900,300));
CHECK(catima::bremsstrahlung(p(1e6),t) == approx(170,20));
CHECK(catima::pair_production(p(7e6),t) == approx(19000,5000));
CHECK(catima::bremsstrahlung(p(7e6),t) == approx(6000,500));
}
TEST_CASE("dEdx for compounds"){
catima::Projectile p{1,1,1,1000};
catima::Material water({
{1.00794,1,2},
{15.9994,8,1}
});
EXPECT( catima::dedx(p(1000), water) == approx(2.23).R(5e-3));
EXPECT( catima::dedx(p(500), water) == approx(2.76).R(5e-3));
EXPECT( catima::dedx(p(9), water) == approx(51.17).R(5e-3));
},
CASE("dEdx from spline vs dEdx"){
CHECK( catima::dedx(p(1000), water) == approx(2.23).R(5e-3));
CHECK( catima::dedx(p(500), water) == approx(2.76).R(5e-3));
CHECK( catima::dedx(p(9), water) == approx(51.17).R(5e-3));
}
TEST_CASE("dEdx from spline vs dEdx"){
catima::Projectile p{238,92,92,1000};
catima::Material graphite({
{12.011,6,1},
});
auto res = catima::calculate(p(1000),graphite);
EXPECT(catima::dedx(p(1000), graphite) == approx(res.dEdxi).R(0.001) );
CHECK(catima::dedx(p(1000), graphite) == approx(res.dEdxi).R(0.001) );
res = catima::calculate(p,graphite,500);
EXPECT(catima::dedx(p(500), graphite) == approx(res.dEdxi).R(0.001) );
CHECK(catima::dedx(p(500), graphite) == approx(res.dEdxi).R(0.001) );
res = catima::calculate(p,graphite,9);
EXPECT(catima::dedx(p(9), graphite) == approx(res.dEdxi).R(0.001) );
},
// CASE("dOmega2dx Ferisov test"){
//},
CASE("Eout test"){
CHECK(catima::dedx(p(9), graphite) == approx(res.dEdxi).R(0.001) );
}
TEST_CASE("Eout test"){
catima::Projectile p{12,6,6,1000};
catima::Material water({
{1.00794,1,2},
@ -167,9 +162,9 @@ const lest::test specification[] =
graphite.thickness(0.5);
auto res = catima::calculate(p,graphite);
EXPECT( res.Eout == approx(997.07,01));
},
CASE("TOF test"){
CHECK( res.Eout == approx(997.07,01));
}
TEST_CASE("TOF test"){
catima::Projectile p{12,6,6,1000};
catima::Material water({
{1.00794,1,2},
@ -186,9 +181,9 @@ const lest::test specification[] =
auto res = catima::calculate(p,water);
dif = res.tof - 0.038;
EXPECT( fabs(dif) < 0.01);
},
CASE("result from stopped material"){
CHECK( fabs(dif) < 0.01);
}
TEST_CASE("result from stopped material"){
catima::Projectile p{12,6,6,1000};
catima::Material water({
{1.00794,1,2},
@ -197,26 +192,26 @@ const lest::test specification[] =
water.density(1.0);
water.thickness(1000.0);
auto res = catima::calculate(p,water);
EXPECT(res.Eout == 0.0);
EXPECT(res.Eloss == 1000*12);
EXPECT(res.sigma_E == 0.0);
EXPECT(res.sigma_a == 0.0);
EXPECT(res.sigma_r > 0.0);
EXPECT(res.dEdxo == 0.0);
EXPECT(res.tof == 0.0);
CHECK(res.Eout == 0.0);
CHECK(res.Eloss == 1000*12);
CHECK(res.sigma_E == 0.0);
CHECK(res.sigma_a == 0.0);
CHECK(res.sigma_r > 0.0);
CHECK(res.dEdxo == 0.0);
CHECK(res.tof == 0.0);
catima::Layers mat;
mat.add(water);
auto res2= catima::calculate(p,mat);
EXPECT(res2.results.size() == 1);
EXPECT(res2.total_result.Eout == res2.results[0].Eout);
EXPECT(res2.total_result.Eout == 0.0);
EXPECT(res2.total_result.Eloss == 1000*12);
EXPECT(res2.total_result.sigma_E == 0.0);
EXPECT(res2.total_result.sigma_a == 0.0);
EXPECT(res2.total_result.tof == 0.0);
},
CASE("constant results from material"){
CHECK(res2.results.size() == 1);
CHECK(res2.total_result.Eout == res2.results[0].Eout);
CHECK(res2.total_result.Eout == 0.0);
CHECK(res2.total_result.Eloss == 1000*12);
CHECK(res2.total_result.sigma_E == 0.0);
CHECK(res2.total_result.sigma_a == 0.0);
CHECK(res2.total_result.tof == 0.0);
}
TEST_CASE("constant results from material"){
catima::Projectile p{12,6,6,1000};
catima::Material water({
{1.00794,1,2},
@ -226,15 +221,15 @@ const lest::test specification[] =
water.thickness(10.0);
auto res = catima::calculate(p,water);
auto res2 = catima::calculate(p,water);
EXPECT(res.Eout == res2.Eout);
EXPECT(res.Eloss == res2.Eloss);
EXPECT(res.sigma_E == res2.sigma_E);
EXPECT(res.sigma_a == res2.sigma_a);
EXPECT(res.sigma_r == res2.sigma_r);
EXPECT(res.dEdxo == res2.dEdxo);
EXPECT(res.tof == res2.tof);
},
CASE("simplified calculation"){
CHECK(res.Eout == res2.Eout);
CHECK(res.Eloss == res2.Eloss);
CHECK(res.sigma_E == res2.sigma_E);
CHECK(res.sigma_a == res2.sigma_a);
CHECK(res.sigma_r == res2.sigma_r);
CHECK(res.dEdxo == res2.dEdxo);
CHECK(res.tof == res2.tof);
}
TEST_CASE("simplified calculation"){
catima::Projectile p{12,6,6,1000};
catima::Material graphite({
{12.011,6,1},
@ -242,28 +237,28 @@ const lest::test specification[] =
graphite.density(2.0).thickness(1.0);
auto res1 = catima::calculate(p,graphite);
auto res2 = catima::calculate(12,6,1000,12.011,6,1.0,2.0);
EXPECT(res1.Eout == res2.Eout);
EXPECT(res1.Eloss == res2.Eloss);
EXPECT(res1.sigma_E == res2.sigma_E);
EXPECT(res1.sigma_a == res2.sigma_a);
EXPECT(res1.sigma_r == res2.sigma_r);
EXPECT(res1.dEdxo == res2.dEdxo);
EXPECT(res1.tof == res2.tof);
CHECK(res1.Eout == res2.Eout);
CHECK(res1.Eloss == res2.Eloss);
CHECK(res1.sigma_E == res2.sigma_E);
CHECK(res1.sigma_a == res2.sigma_a);
CHECK(res1.sigma_r == res2.sigma_r);
CHECK(res1.dEdxo == res2.dEdxo);
CHECK(res1.tof == res2.tof);
auto ra = catima::angular_straggling_from_E(p,res1.Ein,res1.Eout,graphite);
EXPECT(res1.sigma_a == ra);
CHECK(res1.sigma_a == ra);
auto re = catima::energy_straggling_from_E(p,res1.Ein,res1.Eout,graphite);
EXPECT(res1.sigma_E == re);
CHECK(res1.sigma_E == re);
auto eo1 = catima::energy_out(p(1000),graphite);
EXPECT(res1.Eout == eo1);
CHECK(res1.Eout == eo1);
auto de1 = catima::dedx_from_range(p(1000),graphite);
EXPECT(res1.dEdxi == de1);
CHECK(res1.dEdxi == de1);
},
CASE("multilayer basic"){
}
TEST_CASE("multilayer basic"){
catima::Projectile p{12,6,6,1000};
catima::Material water({
{1.00794,1,2},
@ -282,47 +277,47 @@ const lest::test specification[] =
mat.add(graphite);
auto res = catima::calculate(p(1000),mat);
EXPECT(res.total_result.Eout == approx(926.3,0.1));
EXPECT(res.total_result.sigma_a == approx(0.00269).R(0.05));
EXPECT(res.total_result.tof == approx(0.402).R(0.001));
EXPECT(res.total_result.Eloss == approx(884.2,1.0));
//EXPECT(rcompare(res.total_result.sigma_E,0.7067,0.2));
EXPECT(res.results[0].Eout == approx(932.24,0.1));
EXPECT(res.results[0].sigma_a == approx(0.00258).R(0.05));
EXPECT(res.results[0].range == approx(107.163,0.1));
EXPECT(res.results[1].Eout == approx(926.3,0.1));
EXPECT(res.results[1].sigma_a == approx(0.000774).R(0.05));
EXPECT(res.results[1].range == approx(111.3,0.1));
CHECK(res.total_result.Eout == approx(926.3,0.1));
CHECK(res.total_result.sigma_a == approx(0.00269).R(0.05));
CHECK(res.total_result.tof == approx(0.402).R(0.001));
CHECK(res.total_result.Eloss == approx(884.2,1.0));
//CHECK(rcompare(res.total_result.sigma_E,0.7067,0.2));
CHECK(res.results[0].Eout == approx(932.24,0.1));
CHECK(res.results[0].sigma_a == approx(0.00258).R(0.05));
CHECK(res.results[0].range == approx(107.163,0.1));
CHECK(res.results[1].Eout == approx(926.3,0.1));
CHECK(res.results[1].sigma_a == approx(0.000774).R(0.05));
CHECK(res.results[1].range == approx(111.3,0.1));
auto res0 = catima::calculate(p(1000),water);
EXPECT(res0.Eout == res.results[0].Eout);
EXPECT(res0.sigma_a == res.results[0].sigma_a);
EXPECT(res0.sigma_E == res.results[0].sigma_E);
EXPECT(res0.sigma_r == res.results[0].sigma_r);
EXPECT(res0.tof == res.results[0].tof);
CHECK(res0.Eout == res.results[0].Eout);
CHECK(res0.sigma_a == res.results[0].sigma_a);
CHECK(res0.sigma_E == res.results[0].sigma_E);
CHECK(res0.sigma_r == res.results[0].sigma_r);
CHECK(res0.tof == res.results[0].tof);
},
CASE("default material calculations"){
}
TEST_CASE("default material calculations"){
catima::Projectile p{12,6,6,350};
auto air = catima::get_material(catima::material::Air);
air.thickness(0.500);
auto res = catima::calculate(p(350),air);
EXPECT(res.Eout == approx(345.6).epsilon(1.0));
EXPECT(res.sigma_a == approx(0.0013).epsilon(1e-4));
EXPECT(res.sigma_E == approx(0.12).epsilon(1e-3));
EXPECT(res.dEdxi == approx(103.5).epsilon(1e-1));
CHECK(res.Eout == approx(345.6).epsilon(1.0));
CHECK(res.sigma_a == approx(0.0013).epsilon(1e-4));
CHECK(res.sigma_E == approx(0.12).epsilon(1e-3));
CHECK(res.dEdxi == approx(103.5).epsilon(1e-1));
res = catima::calculate(p(150),air);
EXPECT(res.dEdxi == approx(173.6).epsilon(1e0));
CHECK(res.dEdxi == approx(173.6).epsilon(1e0));
res = catima::calculate(p(1000),air);
EXPECT(res.dEdxi == approx(70.69).epsilon(1e-0));
CHECK(res.dEdxi == approx(70.69).epsilon(1e-0));
auto water = catima::get_material(catima::material::Water);
auto res2 = catima::calculate(p(600),water,600);
EXPECT(res2.dEdxi == approx(92.5).epsilon(2));
},
CASE("z_eff"){
CHECK(res2.dEdxi == approx(92.5).epsilon(2));
}
TEST_CASE("z_eff"){
using namespace catima;
Projectile p_u(238,92);
Target t;
@ -330,41 +325,41 @@ const lest::test specification[] =
Config c;
c.z_effective = z_eff_type::pierce_blann;
EXPECT(z_eff_Pierce_Blann(92,beta_from_T(5000.)) == approx(91.8).epsilon(0.2));
EXPECT(z_eff_Pierce_Blann(92,beta_from_T(5000.)) == z_effective(p_u(5000.),t,c));
CHECK(z_eff_Pierce_Blann(92,beta_from_T(5000.)) == approx(91.8).epsilon(0.2));
CHECK(z_eff_Pierce_Blann(92,beta_from_T(5000.)) == z_effective(p_u(5000.),t,c));
EXPECT(z_eff_Winger(92,0.99,6) == approx(91.8).epsilon(0.5));
EXPECT(z_eff_Winger(92,beta_from_T(5000.),13) == approx(91.8).epsilon(0.2));
CHECK(z_eff_Winger(92,0.99,6) == approx(91.8).epsilon(0.5));
CHECK(z_eff_Winger(92,beta_from_T(5000.),13) == approx(91.8).epsilon(0.2));
c.z_effective = z_eff_type::winger;
EXPECT(z_eff_Winger(92,beta_from_T(5000.),13) == z_effective(p_u(5000.),t,c));
CHECK(z_eff_Winger(92,beta_from_T(5000.),13) == z_effective(p_u(5000.),t,c));
EXPECT(z_eff_Schiwietz(92,0.99,6) == approx(91.8).epsilon(0.5));
CHECK(z_eff_Schiwietz(92,0.99,6) == approx(91.8).epsilon(0.5));
c.z_effective = z_eff_type::schiwietz;
EXPECT(z_eff_Schiwietz(92,beta_from_T(5000.),13) == z_effective(p_u(5000.),t,c));
CHECK(z_eff_Schiwietz(92,beta_from_T(5000.),13) == z_effective(p_u(5000.),t,c));
EXPECT(z_eff_Hubert(92,1900,13) == approx(91.88).epsilon(0.1));
CHECK(z_eff_Hubert(92,1900,13) == approx(91.88).epsilon(0.1));
c.z_effective = z_eff_type::hubert;
EXPECT(z_eff_Hubert(92,1900,13) == z_effective(p_u(1900.),t,c));
CHECK(z_eff_Hubert(92,1900,13) == z_effective(p_u(1900.),t,c));
#ifdef GLOBAL
EXPECT(z_eff_global(92,1900,13) == approx(91.88).epsilon(0.05));
CHECK(z_eff_global(92,1900,13) == approx(91.88).epsilon(0.05));
c.z_effective = z_eff_type::global;
EXPECT(z_eff_global(92,1900,13) == z_effective(p_u(1900.),t,c));
EXPECT(z_eff_global(92,1000,13) == approx(91.71).epsilon(0.05));
EXPECT(z_eff_global(92,500,13) == approx(91.22).epsilon(0.1));
EXPECT(z_eff_global(92,100,6) == approx(89.61).epsilon(0.2));
//EXPECT(z_eff_global(92,100,13) == approx(89.42).epsilon(0.1));
//EXPECT(z_eff_global(92,100,29) == approx(88.37).epsilon(0.1));
//EXPECT(z_eff_global(92,50,13) == approx(85.94).epsilon(0.1));
EXPECT(z_eff_global(92,2001,13) == approx(92.0).epsilon(0.01));
EXPECT(z_eff_global(92,2000,13) == approx(92.0).epsilon(0.2));
CHECK(z_eff_global(92,1900,13) == z_effective(p_u(1900.),t,c));
CHECK(z_eff_global(92,1000,13) == approx(91.71).epsilon(0.05));
CHECK(z_eff_global(92,500,13) == approx(91.22).epsilon(0.1));
CHECK(z_eff_global(92,100,6) == approx(89.61).epsilon(0.2));
//CHECK(z_eff_global(92,100,13) == approx(89.42).epsilon(0.1));
//CHECK(z_eff_global(92,100,29) == approx(88.37).epsilon(0.1));
//CHECK(z_eff_global(92,50,13) == approx(85.94).epsilon(0.1));
CHECK(z_eff_global(92,2001,13) == approx(92.0).epsilon(0.01));
CHECK(z_eff_global(92,2000,13) == approx(92.0).epsilon(0.2));
EXPECT(z_eff_atima14(92,1900,13) == approx(91.88).epsilon(0.05));
CHECK(z_eff_atima14(92,1900,13) == approx(91.88).epsilon(0.05));
c.z_effective = z_eff_type::atima14;
EXPECT(z_eff_atima14(92,1900,13) == z_effective(p_u(1900.),t,c));
CHECK(z_eff_atima14(92,1900,13) == z_effective(p_u(1900.),t,c));
#endif
},
CASE("vector_inputs"){
}
TEST_CASE("vector_inputs"){
catima::Projectile p{12,6,6,1000};
catima::Material water({
{1.00794,1,2},
@ -376,31 +371,25 @@ const lest::test specification[] =
graphite.thickness(0.5);
auto res = catima::calculate(p,graphite);
EXPECT( res.Eout == approx(997.07,01));
CHECK( res.Eout == approx(997.07,01));
std::vector<double> energies{100,500,1000};
auto res2 = catima::energy_out(p,energies, graphite);
EXPECT(res2.size()==energies.size());
EXPECT(res2[2] == approx(997.07,01));
EXPECT(res2[0] == approx(catima::energy_out(p(energies[0]),graphite),0.1));
EXPECT(res2[1] == approx(catima::energy_out(p(energies[1]),graphite),0.1));
EXPECT(res2[2] == approx(catima::energy_out(p(energies[2]),graphite),0.1));
CHECK(res2.size()==energies.size());
CHECK(res2[2] == approx(997.07,01));
CHECK(res2[0] == approx(catima::energy_out(p(energies[0]),graphite),0.1));
CHECK(res2[1] == approx(catima::energy_out(p(energies[1]),graphite),0.1));
CHECK(res2[2] == approx(catima::energy_out(p(energies[2]),graphite),0.1));
auto res3 = catima::dedx_from_range(p,energies,graphite);
EXPECT(res3.size()==energies.size());
EXPECT(res3[0] == approx(catima::dedx_from_range(p(energies[0]),graphite),0.1));
EXPECT(res3[1] == approx(catima::dedx_from_range(p(energies[1]),graphite),0.1));
EXPECT(res3[2] == approx(catima::dedx_from_range(p(energies[2]),graphite),0.1));
},
CASE("constants"){
CHECK(res3.size()==energies.size());
CHECK(res3[0] == approx(catima::dedx_from_range(p(energies[0]),graphite),0.1));
CHECK(res3[1] == approx(catima::dedx_from_range(p(energies[1]),graphite),0.1));
CHECK(res3[2] == approx(catima::dedx_from_range(p(energies[2]),graphite),0.1));
}
TEST_CASE("constants"){
using namespace catima;
EXPECT(0.1*hbar*c_light/atomic_mass_unit == approx(0.21183,0.0001));
EXPECT(16.0*dedx_constant*electron_mass*fine_structure/(atomic_mass_unit*3.0*4.0*M_PI) == approx(5.21721169334564e-7).R(1e-3));
CHECK(0.1*hbar*c_light/atomic_mass_unit == approx(0.21183,0.0001));
CHECK(16.0*dedx_constant*electron_mass*fine_structure/(atomic_mass_unit*3.0*4.0*M_PI) == approx(5.21721169334564e-7).R(1e-3));
}
};
int main( int argc, char * argv[] )
{
return lest::run( specification, argc, argv );
}

View File

@ -1,25 +1,14 @@
#include "lest.hpp"
#define DOCTEST_CONFIG_IMPLEMENT_WITH_MAIN
#define DOCTEST_CONFIG_SUPER_FAST_ASSERTS
#include "doctest.h"
#include <math.h>
#include "testutils.h"
#include "catima/catima.h"
//#include "nucdata.h"
using namespace std;
using lest::approx;
bool rcompare(double a, double b,double eps){
if(fabs((a-b)/fabs(b))<eps){
return true;
}
else{
std::cout<<"\033[1;31m"<<a<<" == "<<b<<"\033[0m"<<std::endl;
return false;
}
}
const lest::test specification[] =
{
CASE("LS generated is equal to calculated"){
TEST_CASE("LS generated is equal to calculated"){
catima::Projectile p;
double a,b;
@ -29,7 +18,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard(p);
b = catima::precalculated_lindhard(p);
EXPECT(a==approx(b).epsilon(0.001));
CHECK(a==approx(b).epsilon(0.001));
}
p.A = 220;
@ -38,7 +27,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard(p);
b = catima::precalculated_lindhard(p);
EXPECT(a==approx(b).epsilon(0.01));
CHECK(a==approx(b).epsilon(0.01));
}
p.A = 250;
@ -47,7 +36,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard(p);
b = catima::precalculated_lindhard(p);
EXPECT(a==approx(b).epsilon(0.01));
CHECK(a==approx(b).epsilon(0.01));
}
@ -57,7 +46,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard(p);
b = catima::precalculated_lindhard(p);
EXPECT(a==approx(b).epsilon(0.01));
CHECK(a==approx(b).epsilon(0.01));
}
p.A = 100;
@ -66,12 +55,12 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard(p);
b = catima::precalculated_lindhard(p);
EXPECT(a==approx(b).epsilon(0.01));
CHECK(a==approx(b).epsilon(0.01));
}
},
}
CASE("LS X generated is equal to calculated"){
TEST_CASE("LS X generated is equal to calculated"){
catima::Projectile p;
double a,b;
@ -81,7 +70,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard_X(p);
b = catima::precalculated_lindhard_X(p);
EXPECT(a==approx(b).epsilon(0.001));
CHECK(a==approx(b).epsilon(0.001));
}
p.A = 220;
@ -90,7 +79,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard_X(p);
b = catima::precalculated_lindhard_X(p);
EXPECT(a==approx(b).epsilon(0.02));
CHECK(a==approx(b).epsilon(0.02));
}
p.A = 250;
@ -99,7 +88,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard_X(p);
b = catima::precalculated_lindhard_X(p);
EXPECT(a==approx(b).epsilon(0.03));
CHECK(a==approx(b).epsilon(0.03));
}
@ -109,7 +98,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard_X(p);
b = catima::precalculated_lindhard_X(p);
EXPECT(a==approx(b).epsilon(0.03));
CHECK(a==approx(b).epsilon(0.03));
}
p.A = 100;
@ -118,15 +107,7 @@ const lest::test specification[] =
p.T = e;
a = catima::bethek_lindhard_X(p);
b = catima::precalculated_lindhard_X(p);
EXPECT(a==approx(b).epsilon(0.03));
CHECK(a==approx(b).epsilon(0.03));
}
}
};
int main( int argc, char * argv[] )
{
return lest::run( specification, argc, argv );
}

View File

@ -1,15 +1,14 @@
#include "lest.hpp"
#define DOCTEST_CONFIG_IMPLEMENT_WITH_MAIN
#define DOCTEST_CONFIG_SUPER_FAST_ASSERTS
#include "doctest.h"
#include <math.h>
#include "catima/catima.h"
#include "catima/reactions.h"
#include "testutils.h"
using namespace std;
using catima::approx;
using catima::reaction_rate;
using catima::nonreaction_rate;
const lest::test specification[] =
{
CASE("reaction"){
TEST_CASE("reaction"){
catima::Projectile proj{12,6,6,870};
catima::Projectile proj2{238,92,92,500};
auto c = catima::get_material(6);
@ -19,7 +18,7 @@ const lest::test specification[] =
double r,r2;
r= catima::nonreaction_rate(proj, c);
EXPECT(r == approx(0.92,0.02));
CHECK(r == approx(0.92,0.02));
catima::Layers l;
l.add(c);
@ -27,29 +26,29 @@ const lest::test specification[] =
l.add(c);
auto res = catima::calculate(proj,l);
EXPECT(res.total_result.sp == approx(r*r*r,0.02));
CHECK(res.total_result.sp == approx(r*r*r,0.02));
c.thickness(6.0);
r= catima::nonreaction_rate(proj, c);
EXPECT(res.total_result.sp == approx(r,0.001));
CHECK(res.total_result.sp == approx(r,0.001));
c.thickness(0.0);
r= catima::nonreaction_rate(proj, c);
EXPECT(r == 1.0);
CHECK(r == 1.0);
proj.T = 0.0;
c.thickness(6);
r= catima::nonreaction_rate(proj, c);
EXPECT(r == -1.0);
CHECK(r == -1.0);
proj.T=870;
water.thickness(1);
r= catima::nonreaction_rate(proj2, water);
r2= catima::nonreaction_rate(proj, water);
EXPECT( (r > 0 && r<1.0) );
EXPECT( (r2 > 0 && r2<1.0) );
EXPECT( r2>r );
},
CASE("production"){
CHECK( (r > 0 && r<1.0) );
CHECK( (r2 > 0 && r2<1.0) );
CHECK( r2>r );
}
TEST_CASE("production"){
catima::Projectile proj{12,6,6,870};
auto c = catima::get_material(6);
c.thickness(0.1);
@ -61,23 +60,16 @@ const lest::test specification[] =
r = 0.0001*cs*c.number_density_cm2();
r2 = catima::production_rate(cs,rcsi,rcso, c);
EXPECT(r==approx(r2).R(0.01));
CHECK(r==approx(r2).R(0.01));
r2 = catima::production_rate(cs,2,1, c);
EXPECT(r==approx(r2).R(0.001));
CHECK(r==approx(r2).R(0.001));
r = catima::production_rate(cs,870,870, c);
r2 = catima::production_rate(cs,870,860, c);
EXPECT(r==approx(r2).R(0.001));
CHECK(r==approx(r2).R(0.001));
c.thickness(2.0);
r = catima::production_rate(cs,870,870, c);
r2 = catima::production_rate(cs,870,860, c);
EXPECT(r==approx(r2).R(0.001));
CHECK(r==approx(r2).R(0.001));
}
};
int main( int argc, char * argv[] )
{
return lest::run( specification, argc, argv );
}

View File

@ -1,4 +1,6 @@
#include "lest.hpp"
#define DOCTEST_CONFIG_IMPLEMENT_WITH_MAIN
#define DOCTEST_CONFIG_SUPER_FAST_ASSERTS
#include "doctest.h"
#include <math.h>
#include "testutils.h"
using namespace std;
@ -6,11 +8,7 @@ using namespace std;
#include "catima/catima.h"
#include "catima/storage.h"
using catima::approx;
const lest::test specification[] =
{
CASE("datapoints equal operator"){
TEST_CASE("datapoints equal operator"){
catima::Projectile p{12,6,6,1000};
catima::Material water({
{1,1,2},
@ -28,18 +26,18 @@ const lest::test specification[] =
catima::DataPoint d;
catima::DataPoint e(p,water,c2);
d = c;
EXPECT(a == b);
EXPECT(!(a==c));
EXPECT(d == c);
EXPECT(!(d==b));
EXPECT(!(e==a));
CHECK(a == b);
CHECK(!(a==c));
CHECK(d == c);
CHECK(!(d==b));
CHECK(!(e==a));
d = a;
EXPECT(!(d==c));
EXPECT(d==b);
CHECK(!(d==c));
CHECK(d==b);
},
}
CASE("storage add"){
TEST_CASE("storage add"){
catima::Projectile p{12,6,6,1000};
catima::Material water({
{1,1,2},
@ -51,42 +49,42 @@ const lest::test specification[] =
});
catima::_storage.Reset();
EXPECT(catima::_storage.get_index()==0);
CHECK(catima::_storage.get_index()==0);
catima::_storage.Add(p,water);
auto& dp = catima::_storage.Get(0);
EXPECT(catima::_storage.get_index()==1);
EXPECT(dp.p.A==12);
EXPECT(dp.m.ncomponents()==2);
CHECK(catima::_storage.get_index()==1);
CHECK(dp.p.A==12);
CHECK(dp.m.ncomponents()==2);
catima::_storage.Add(p,water);
auto& dp2 = catima::_storage.Get(1);
EXPECT(catima::_storage.get_index()==1);
EXPECT(dp2.p.A==0);
EXPECT(dp2.m.ncomponents()==0);
CHECK(catima::_storage.get_index()==1);
CHECK(dp2.p.A==0);
CHECK(dp2.m.ncomponents()==0);
catima::_storage.Add(p,graphite);
auto& dp3 = catima::_storage.Get(1);
EXPECT(catima::_storage.get_index()==2);
EXPECT(dp3.p.A==12);
EXPECT(dp3.m.ncomponents()==1);
CHECK(catima::_storage.get_index()==2);
CHECK(dp3.p.A==12);
CHECK(dp3.m.ncomponents()==1);
catima::_storage.Add(p,graphite);
EXPECT(catima::_storage.get_index()==2);
CHECK(catima::_storage.get_index()==2);
catima::Config c1;
c1.z_effective = catima::z_eff_type::global;
catima::_storage.Add(p,graphite, c1);
EXPECT(catima::_storage.get_index()==3);
CHECK(catima::_storage.get_index()==3);
catima::_storage.Add(p,graphite);
EXPECT(catima::_storage.get_index()==3);
CHECK(catima::_storage.get_index()==3);
c1.z_effective = catima::z_eff_type::hubert;
catima::_storage.Add(p,graphite ,c1);
EXPECT(catima::_storage.get_index()==4);
CHECK(catima::_storage.get_index()==4);
},
CASE("test maximum storage"){
}
TEST_CASE("test maximum storage"){
auto maxdata = catima::max_storage_data;
catima::Projectile p{12,6,6,1000};
catima::Material water({
@ -98,53 +96,45 @@ const lest::test specification[] =
{12,6,1}
});
catima::_storage.Reset();
EXPECT(catima::_storage.get_index()==0);
CHECK(catima::_storage.get_index()==0);
for(int i=1;i<maxdata+1;i++){
catima::Projectile p1{2.0*i,(double)i,(double)i,1000};
catima::_storage.Add(p1,graphite);
EXPECT(catima::_storage.get_index()==i);
EXPECT(catima::_storage.GetN()==maxdata);
CHECK(catima::_storage.get_index()==i);
CHECK(catima::_storage.GetN()==maxdata);
}
EXPECT(catima::_storage.get_index()==maxdata);
CHECK(catima::_storage.get_index()==maxdata);
for(int i=1;i<maxdata-1;i++){
catima::Projectile p1{2.0*i,(double)i,(double)i,1000};
catima::_storage.Add(p1,water);
EXPECT(catima::_storage.get_index()==i);
EXPECT(catima::_storage.GetN()==maxdata);
CHECK(catima::_storage.get_index()==i);
CHECK(catima::_storage.GetN()==maxdata);
}
},
CASE("energy table"){
}
TEST_CASE("energy table"){
double step = (catima::logEmax - catima::logEmin)/(catima::max_datapoints-1);
EXPECT(catima::energy_table.step==step);
EXPECT(catima::energy_table.values[0]==exp(M_LN10*(catima::logEmin)));
EXPECT(catima::energy_table.values[1]==exp(M_LN10*(catima::logEmin+step)));
EXPECT(catima::energy_table.values[2]==exp(M_LN10*(catima::logEmin+2.0*step)));
EXPECT(catima::energy_table.values[3]==exp(M_LN10*(catima::logEmin+3.0*step)));
EXPECT(catima::energy_table.values[4]==exp(M_LN10*(catima::logEmin+4.0*step)));
EXPECT(catima::energy_table.values[5]==exp(M_LN10*(catima::logEmin+5.0*step)));
EXPECT(catima::energy_table.values[catima::max_datapoints-1]==approx(exp(M_LN10*(catima::logEmax))).epsilon(1e-6));
},
CASE("indexing"){
CHECK(catima::energy_table.step==step);
CHECK(catima::energy_table.values[0]==exp(M_LN10*(catima::logEmin)));
CHECK(catima::energy_table.values[1]==exp(M_LN10*(catima::logEmin+step)));
CHECK(catima::energy_table.values[2]==exp(M_LN10*(catima::logEmin+2.0*step)));
CHECK(catima::energy_table.values[3]==exp(M_LN10*(catima::logEmin+3.0*step)));
CHECK(catima::energy_table.values[4]==exp(M_LN10*(catima::logEmin+4.0*step)));
CHECK(catima::energy_table.values[5]==exp(M_LN10*(catima::logEmin+5.0*step)));
CHECK(catima::energy_table.values[catima::max_datapoints-1]==approx(exp(M_LN10*(catima::logEmax))).epsilon(1e-6));
}
TEST_CASE("indexing"){
double val, dif;
EXPECT(EnergyTable_index(catima::energy_table, 0.0)==-1);
CHECK(EnergyTable_index(catima::energy_table, 0.0)==-1);
for(int i=0;i<catima::max_datapoints-1;i++){
val = catima::energy_table.values[i];
dif = catima::energy_table.values[i+1] - val;
EXPECT(EnergyTable_index(catima::energy_table, val)==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);
CHECK(EnergyTable_index(catima::energy_table, val)==i);
CHECK(EnergyTable_index(catima::energy_table, val+0.5*dif)==i);
CHECK(catima::energy_table.index(val)==i);
CHECK(catima::energy_table.index(val+0.5*dif)==i);
}
}
};
int main( int argc, char * argv[] )
{
return lest::run( specification, argc, argv );
}
}

View File

@ -1,9 +1,9 @@
#include "lest.hpp"
#define DOCTEST_CONFIG_IMPLEMENT_WITH_MAIN
#define DOCTEST_CONFIG_SUPER_FAST_ASSERTS
#include "doctest.h"
#include "testutils.h"
#include <math.h>
using namespace std;
using catima::approx;
#include "catima/catima.h"
#include "catima/material_database.h"
@ -18,10 +18,7 @@ bool rcompare(double a, double b,double eps){
}
const lest::test specification[] =
{
CASE("atima material basic tests"){
SETUP(""){
TEST_CASE("atima material basic tests"){
catima::Material water({
{1,1,2},
{16,8,1}
@ -35,63 +32,54 @@ const lest::test specification[] =
graphite.add_element(12,6,1);
graphite.density(1.8);
SECTION("component number check"){
EXPECT(graphite.ncomponents()==1);
EXPECT(water.ncomponents()==2);
graphite.add_element(18,40,1);
EXPECT(graphite.ncomponents()==2);
EXPECT(water.M()==18);
}
SECTION("Molar mass check"){
EXPECT(graphite.M()==12);
EXPECT(water.ncomponents()==2);
EXPECT(water.M()==18);
EXPECT(water2.ncomponents()==2);
EXPECT(water2.M()==18);
}
SECTION("equal operator check"){
EXPECT(water==water2);
EXPECT(!(water==graphite));
}
SECTION("default ionisation potential"){
EXPECT(graphite.I()==0.0);
}
SECTION("length"){
water.density(1.0);
water.thickness(1.0);
EXPECT(water.thickness()==approx(1.0,0.0001));
water.thickness_cm(1.0);
EXPECT(water.thickness()==approx(1.0,0.0001));
water.thickness_cm(2.0);
EXPECT(water.thickness()==approx(2.0,0.0001));
}
}
},
CASE("Material automatic atomic weight"){
CHECK(graphite.ncomponents()==1);
CHECK(water.ncomponents()==2);
CHECK(graphite.M()==12);
CHECK(water.ncomponents()==2);
CHECK(water.M()==18);
CHECK(water2.ncomponents()==2);
CHECK(water2.M()==18);
CHECK(graphite.I()==0.0);
graphite.add_element(18,40,1);
CHECK(graphite.ncomponents()==2);
CHECK_FALSE(graphite.M()==approx(12,0.1));
CHECK(water==water2);
CHECK(!(water==graphite));
water.density(1.0);
water.thickness(1.0);
CHECK(water.thickness()==approx(1.0,0.0001));
water.thickness_cm(1.0);
CHECK(water.thickness()==approx(1.0,0.0001));
water.thickness_cm(2.0);
CHECK(water.thickness()==approx(2.0,0.0001));
}
TEST_CASE("Material automatic atomic weight"){
catima::Material water({{0,1,2},{0,8,1}});
catima::Material graphite(0,6);
EXPECT(water.get_element(0).A == 1.00794);
EXPECT(water.get_element(1).A == 15.9994);
EXPECT(graphite.get_element(0).A == 12.0107);
EXPECT(water.M()>16);
EXPECT(graphite.M()>12);
},
CASE("default materials"){
CHECK(water.get_element(0).A == 1.00794);
CHECK(water.get_element(1).A == 15.9994);
CHECK(graphite.get_element(0).A == 12.0107);
CHECK(water.M()>16);
CHECK(graphite.M()>12);
}
TEST_CASE("default materials"){
catima::Material m = catima::get_material(6);
EXPECT(m.get_element(0).A == 12.0107);
EXPECT(m.get_element(0).Z == 6);
EXPECT(m.density() == 2.0);
EXPECT(m.M() == 12.0107);
CHECK(m.get_element(0).A == 12.0107);
CHECK(m.get_element(0).Z == 6);
CHECK(m.density() == 2.0);
CHECK(m.M() == 12.0107);
m = catima::get_material(catima::material::Water);
EXPECT(m.get_element(0).A == 1.00794);
EXPECT(m.get_element(0).Z == 1);
EXPECT(m.get_element(1).A == 15.9994);
EXPECT(m.get_element(1).Z == 8);
EXPECT(m.density() == 1.0);
},
CASE("Layers"){
CHECK(m.get_element(0).A == 1.00794);
CHECK(m.get_element(0).Z == 1);
CHECK(m.get_element(1).A == 15.9994);
CHECK(m.get_element(1).Z == 8);
CHECK(m.density() == 1.0);
}
TEST_CASE("Layers"){
catima::Material water2;
water2.add_element(1,1,2);
water2.add_element(16,8,1);
@ -102,83 +90,83 @@ const lest::test specification[] =
graphite.density(1.8).thickness(1.0);
catima::Layers detector1;
EXPECT(detector1.num() == 0);
CHECK(detector1.num() == 0);
detector1.add(graphite);
EXPECT(detector1.num() == 1);
CHECK(detector1.num() == 1);
detector1.add(water2);
detector1.add(graphite);
EXPECT(detector1.num() == 3);
CHECK(detector1.num() == 3);
// check correct density and thickness
EXPECT(detector1[0].density()==1.8);
EXPECT(detector1[0].thickness()==1.0);
EXPECT(detector1[1].density()==1.0);
EXPECT(detector1[1].thickness()==2.0);
EXPECT(detector1[0].get_element(0).Z == 6);
EXPECT(detector1[0].get_element(0).A == 12);
EXPECT(detector1[1].get_element(0).A == 1);
EXPECT(detector1[1].get_element(0).Z == 1);
EXPECT(detector1[2].density() == 1.8);
CHECK(detector1[0].density()==1.8);
CHECK(detector1[0].thickness()==1.0);
CHECK(detector1[1].density()==1.0);
CHECK(detector1[1].thickness()==2.0);
CHECK(detector1[0].get_element(0).Z == 6);
CHECK(detector1[0].get_element(0).A == 12);
CHECK(detector1[1].get_element(0).A == 1);
CHECK(detector1[1].get_element(0).Z == 1);
CHECK(detector1[2].density() == 1.8);
detector1[1].density(1.2);
EXPECT(detector1[1].density()==1.2);
CHECK(detector1[1].density()==1.2);
catima::Layers detector2;
detector2 = detector1;
EXPECT(detector2.num() == 3);
CHECK(detector2.num() == 3);
detector2.add(water2);
detector2.add(graphite);
EXPECT(detector2.num() == 5);
CHECK(detector2.num() == 5);
catima::Layers focal_material = detector1 + detector2;
EXPECT(focal_material.num() == 8);
CHECK(focal_material.num() == 8);
},
CASE("basic projectile tests"){
}
TEST_CASE("basic projectile tests"){
catima::Projectile p{12,6,6,1000};
EXPECT(p.A==12);
EXPECT(p.Z==6);
EXPECT(p.Q==6);
EXPECT(p.T==1000);
CHECK(p.A==12);
CHECK(p.Z==6);
CHECK(p.Q==6);
CHECK(p.T==1000);
catima::Projectile p2(12,6);
EXPECT(p.A==12);
EXPECT(p2.Z==6);
EXPECT(p2.Q==6);
EXPECT(p2.T==0);
CHECK(p.A==12);
CHECK(p2.Z==6);
CHECK(p2.Q==6);
CHECK(p2.T==0);
p2(1000);
EXPECT(p2.T==1000);
CHECK(p2.T==1000);
p2(1000)(500);
EXPECT(p2.T==500);
CHECK(p2.T==500);
catima::Projectile p3(12,6,5);
EXPECT(p==p2);
EXPECT( !(p==p3));
},
CASE("basic config test"){
CHECK(p==p2);
CHECK( !(p==p3));
}
TEST_CASE("basic config test"){
catima::Config c1;
catima::Config c2;
catima::Config c3{catima::z_eff_type::none};
catima::Config c4;
EXPECT(c1.z_effective == catima::default_config.z_effective);
CHECK(c1.z_effective == catima::default_config.z_effective);
EXPECT(c1==c2);
EXPECT( !(c1==c3));
EXPECT(c1==c4);
CHECK(c1==c2);
CHECK( !(c1==c3));
CHECK(c1==c4);
c4.z_effective = catima::z_eff_type::global;
EXPECT(!(c1==c4));
CHECK(!(c1==c4));
auto c5 = c4;
EXPECT(c4==c5);
CHECK(c4==c5);
c4.z_effective = catima::z_eff_type::hubert;
EXPECT(!(c4==c5) );
EXPECT(!(c4==c1));
CHECK(!(c4==c5) );
CHECK(!(c4==c1));
c4.z_effective = catima::default_config.z_effective;
EXPECT(!(c5==c4));
EXPECT((c1==c4));
},
CASE("constructors test"){
CHECK(!(c5==c4));
CHECK((c1==c4));
}
TEST_CASE("constructors test"){
catima::Material mat2(12,6,2.5,0.1);
catima::Material mat3(12.01,6);
catima::Material mat4({{12.01,6,1.0}});
@ -186,26 +174,26 @@ const lest::test specification[] =
{12.01, 6, 1},
{16.00, 8, 2}
});
EXPECT(mat2.ncomponents()==1);
EXPECT(mat3.ncomponents()==1);
EXPECT(mat3.get_element(0).A==12.01);
EXPECT(mat4.ncomponents()==1);
EXPECT(mat4.get_element(0).A==12.01);
EXPECT(mat5.ncomponents()==2);
EXPECT(mat5.get_element(0).A==12.01);
EXPECT(mat5.get_element(0).Z==6);
EXPECT(mat5.get_element(1).A==16.0);
EXPECT(mat5.get_element(1).stn==2);
CHECK(mat2.ncomponents()==1);
CHECK(mat3.ncomponents()==1);
CHECK(mat3.get_element(0).A==12.01);
CHECK(mat4.ncomponents()==1);
CHECK(mat4.get_element(0).A==12.01);
CHECK(mat5.ncomponents()==2);
CHECK(mat5.get_element(0).A==12.01);
CHECK(mat5.get_element(0).Z==6);
CHECK(mat5.get_element(1).A==16.0);
CHECK(mat5.get_element(1).stn==2);
catima::Material mat6;
mat6 = mat5;
EXPECT(mat5==mat6);
EXPECT(mat5.density() == mat6.density());
EXPECT(mat5.ncomponents()==mat6.ncomponents());
EXPECT(mat5.get_element(0).A==mat6.get_element(0).A);
EXPECT(mat5.get_element(1).A==mat6.get_element(1).A);
CHECK(mat5==mat6);
CHECK(mat5.density() == mat6.density());
CHECK(mat5.ncomponents()==mat6.ncomponents());
CHECK(mat5.get_element(0).A==mat6.get_element(0).A);
CHECK(mat5.get_element(1).A==mat6.get_element(1).A);
mat5.add_element(12,6,1);
EXPECT(mat5.ncomponents()==mat6.ncomponents()+1);
CHECK(mat5.ncomponents()==mat6.ncomponents()+1);
// constructor with custom Ipot
catima::Material water1({
@ -216,16 +204,16 @@ const lest::test specification[] =
{1,1,2},
{16,8,1}
},1.0, 78.0);
EXPECT(water1.ncomponents()==2);
EXPECT(water2.ncomponents()==2);
EXPECT(water1.density()==1.0);
EXPECT(water2.density()==1.0);
EXPECT(water1.I()==0.0);
EXPECT(water2.I()==78.0);
EXPECT_NOT(water1==water2);
CHECK(water1.ncomponents()==2);
CHECK(water2.ncomponents()==2);
CHECK(water1.density()==1.0);
CHECK(water2.density()==1.0);
CHECK(water1.I()==0.0);
CHECK(water2.I()==78.0);
CHECK_FALSE(water1==water2);
},
CASE("fraction vs stn init"){
}
TEST_CASE("fraction vs stn init"){
catima::Projectile p{12,6};
catima::Material water1({
{0, 1, 0.111894},
@ -239,12 +227,12 @@ const lest::test specification[] =
water2.thickness(1.0);
auto res1 = catima::calculate(p(600),water1);
auto res2 = catima::calculate(p(600),water2);
EXPECT(res1.dEdxi == approx(res2.dEdxi,0.001));
EXPECT(res1.range == approx(res2.range).R(1e-2));
EXPECT(res1.sigma_a == approx(res2.sigma_a).R(1e-2));
EXPECT(res1.sigma_r == approx(res2.sigma_r).R(1e-2));
},
CASE("fraction calculation"){
CHECK(res1.dEdxi == approx(res2.dEdxi,0.001));
CHECK(res1.range == approx(res2.range).R(1e-2));
CHECK(res1.sigma_a == approx(res2.sigma_a).R(1e-2));
CHECK(res1.sigma_r == approx(res2.sigma_r).R(1e-2));
}
TEST_CASE("fraction calculation"){
catima::Material water1({
{0, 1, 0.111898},
{0, 8, 0.888102}
@ -256,31 +244,31 @@ const lest::test specification[] =
auto air = catima::get_material(catima::material::Air);
EXPECT(water1.weight_fraction(0)==0.111898);
EXPECT(water2.weight_fraction(0)==approx(water1.weight_fraction(0)).R(1e-5));
EXPECT(water1.weight_fraction(1)==0.888102);
EXPECT(water2.weight_fraction(1)==approx(water1.weight_fraction(1)).R(1e-5));
EXPECT(water2.M()==approx(18).epsilon(0.1));
CHECK(water1.weight_fraction(0)==0.111898);
CHECK(water2.weight_fraction(0)==approx(water1.weight_fraction(0)).R(1e-5));
CHECK(water1.weight_fraction(1)==0.888102);
CHECK(water2.weight_fraction(1)==approx(water1.weight_fraction(1)).R(1e-5));
CHECK(water2.M()==approx(18).epsilon(0.1));
EXPECT(water1.M()==approx(6.0,0.1));
EXPECT(water2.M()==approx(18,0.1));
CHECK(water1.M()==approx(6.0,0.1));
CHECK(water2.M()==approx(18,0.1));
EXPECT(water1.molar_fraction(0)==approx(2.0/3.0).R(1e-5));
EXPECT(water2.molar_fraction(0)==approx(2.0).R(1e-5));
EXPECT(water1.molar_fraction(1)==approx(1.0/3.0).R(1e-5));
EXPECT(water2.molar_fraction(1)==approx(1.0).R(1e-5));
EXPECT(water1.molar_fraction(0)/water1.molar_fraction(1)==approx(2.0).R(1e-5));
EXPECT(water2.molar_fraction(0)/water2.molar_fraction(1)==approx(2.0).R(1e-5));
CHECK(water1.molar_fraction(0)==approx(2.0/3.0).R(1e-5));
CHECK(water2.molar_fraction(0)==approx(2.0).R(1e-5));
CHECK(water1.molar_fraction(1)==approx(1.0/3.0).R(1e-5));
CHECK(water2.molar_fraction(1)==approx(1.0).R(1e-5));
CHECK(water1.molar_fraction(0)/water1.molar_fraction(1)==approx(2.0).R(1e-5));
CHECK(water2.molar_fraction(0)/water2.molar_fraction(1)==approx(2.0).R(1e-5));
catima::Material mat({12.0,6,1});
EXPECT(mat.M()==approx(12.0,0.001));
EXPECT(mat.weight_fraction(0)==approx(1.0).R(1e-6));
CHECK(mat.M()==approx(12.0,0.001));
CHECK(mat.weight_fraction(0)==approx(1.0).R(1e-6));
//EXPECT(air.M() == approx(28.97,0.1));
//CHECK(air.M() == approx(28.97,0.1));
},
CASE("number density"){
}
TEST_CASE("number density"){
catima::Material c({12.0,6,1});
auto water = catima::get_material(catima::material::Water);
auto air = catima::get_material(catima::material::Air);
@ -289,36 +277,28 @@ const lest::test specification[] =
air.density(1.2041e-3);
c.thickness_cm(1.0);
EXPECT(c.number_density()==approx(1.7662,0.01));
EXPECT(c.number_density_cm2()==approx(1.7662,0.01));
EXPECT(c.number_density(0)==approx(1.7662,0.01));
EXPECT(c.number_density_cm2(0)==approx(1.7662,0.01));
EXPECT(c.number_density(1)==0.0);
EXPECT(c.number_density_cm2(1)==0.0);
CHECK(c.number_density()==approx(1.7662,0.01));
CHECK(c.number_density_cm2()==approx(1.7662,0.01));
CHECK(c.number_density(0)==approx(1.7662,0.01));
CHECK(c.number_density_cm2(0)==approx(1.7662,0.01));
CHECK(c.number_density(1)==0.0);
CHECK(c.number_density_cm2(1)==0.0);
c.thickness_cm(2.0);
EXPECT(c.number_density()==approx(1.7662,0.01));
EXPECT(c.number_density_cm2()==approx(2.0*1.7662,0.01));
CHECK(c.number_density()==approx(1.7662,0.01));
CHECK(c.number_density_cm2()==approx(2.0*1.7662,0.01));
water.thickness_cm(1.0);
EXPECT(water.number_density()==approx(0.3336,0.001));
EXPECT(water.number_density_cm2()==approx(0.3336,0.001));
EXPECT(water.number_density(0)==approx(2*0.3336,0.001));
EXPECT(water.number_density_cm2(0)==approx(2*0.3336,0.001));
EXPECT(water.number_density(1)==approx(0.3336,0.001));
EXPECT(water.number_density_cm2(1)==approx(0.3336,0.001));
CHECK(water.number_density()==approx(0.3336,0.001));
CHECK(water.number_density_cm2()==approx(0.3336,0.001));
CHECK(water.number_density(0)==approx(2*0.3336,0.001));
CHECK(water.number_density_cm2(0)==approx(2*0.3336,0.001));
CHECK(water.number_density(1)==approx(0.3336,0.001));
CHECK(water.number_density_cm2(1)==approx(0.3336,0.001));
water.thickness_cm(3.0);
EXPECT(water.number_density_cm2()==approx(3.0*0.3336,0.001));
CHECK(water.number_density_cm2()==approx(3.0*0.3336,0.001));
air.thickness_cm(1.0);
EXPECT(air.number_density(0)==approx(air.molar_fraction(0)*2*0.0002504,0.00001));
EXPECT(air.number_density(1)==approx(air.molar_fraction(1)*2*0.0002504,0.00001));
EXPECT(air.number_density(2)==approx(air.molar_fraction(2)*1*0.0002504,0.00001));
CHECK(air.number_density(0)==approx(air.molar_fraction(0)*2*0.0002504,0.00001));
CHECK(air.number_density(1)==approx(air.molar_fraction(1)*2*0.0002504,0.00001));
CHECK(air.number_density(2)==approx(air.molar_fraction(2)*1*0.0002504,0.00001));
}
};
int main( int argc, char * argv[] )
{
return lest::run( specification, argc, argv );
}

View File

@ -1,5 +1,6 @@
namespace catima{
#include <cmath>
#include <stdio.h>
#include <iostream>
class approx
{
public:
@ -20,7 +21,6 @@ public:
approx & epsilon( double epsilon ) { epsilon_ = epsilon; return *this; }
approx & R( double relative ) { epsilon_ = relative*magnitude_; return *this; }
friend bool operator == ( double lhs, approx const & rhs )
{
return std::abs( lhs - rhs.magnitude_) < rhs.epsilon_;
@ -41,31 +41,5 @@ public:
};
std::ostream & operator<<(std::ostream &os, approx const &a){
using lest::to_string;
return os<<to_string(a.magnitude_)<<" +- "<<a.epsilon_;
}
bool rdiff(double a, double b,double eps){
if(fabs((a-b)/fabs(b))<eps){
return true;
}
else{
std::cout<<"\033[1;31m"<<a<<" == "<<b<<"\033[0m"<<std::endl;
return false;
}
}
bool diff(double a, double b,double eps){
if(fabs((a-b))<eps){
return true;
}
else{
std::cout<<"\033[1;31m"<<a<<" == "<<b<<"\033[0m"<<std::endl;
return false;
}
}
return os<<a.magnitude_<<" +- "<<a.epsilon_;
}