From 2fdf06a45b4a23d56833a5440f1437bd77bb9d5a Mon Sep 17 00:00:00 2001 From: hrocho Date: Fri, 19 Oct 2018 02:51:01 +0200 Subject: [PATCH] partial update --- .gitignore | 1 + CMakeLists.txt | 23 ++++----- integrator.cpp | 98 ------------------------------------- integrator.h | 100 +++++++++++++++++++++++++++++++++----- tests/test.py | 17 +++++++ tests/test_dedx_range.cpp | 6 +-- 6 files changed, 122 insertions(+), 123 deletions(-) diff --git a/.gitignore b/.gitignore index 8d18470..63a6f81 100644 --- a/.gitignore +++ b/.gitignore @@ -35,4 +35,5 @@ ana/* benchmark/* build/* +global/* .vscode/* diff --git a/CMakeLists.txt b/CMakeLists.txt index 06b831d..86f506c 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -15,19 +15,24 @@ option(REACTIONS "enable/disable nuclear reaction rate" ON) option(APPS "build catima applications" ON) ######## build type ############ -set(CMAKE_BUILD_TYPE Release) -#set(CMAKE_BUILD_TYPE Debug) -#set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} -g") +if(NOT CMAKE_BUILD_TYPE STREQUAL "Debug") + set(CMAKE_BUILD_TYPE "Release") + MESSAGE(STATUS "Build type Release") +else() + set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} -g") + if (${CMAKE_CXX_COMPILER_ID} STREQUAL "Clang" OR ${CMAKE_CXX_COMPILER_ID} STREQUAL "GNU") + set(CMAKE_CXX_FLAGS_DEBUG "-Wall -Wextra -Wfatal-errors -g") + endif() + MESSAGE(STATUS "Build type Debug") +endif() +MESSAGE(STATUS "Build type: " ${CMAKE_BUILD_TYPE}) ################################ ######### compiler flags ########### set(CMAKE_CXX_EXTENSIONS OFF) set(CMAKE_CXX_STANDARD 14) set(CMAKE_CXX_STANDARD_REQUIRED ON) -if ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "GNU" OR - "${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang") - set(CMAKE_CXX_FLAGS "-Wall -Wextra -Werror -Wno-sign-compare -Wno-unused-parameter") -endif() +MESSAGE(STATUS "install prefix: " ${CMAKE_INSTALL_PREFIX}) if(APPLE) set(RPATH_VARIABLE "DYLD_LIBRARY_PATH") @@ -35,8 +40,6 @@ else() set(RPATH_VARIABLE "LD_LIBRARY_PATH") endif() -MESSAGE(STATUS "install prefix: " ${CMAKE_INSTALL_PREFIX}) - ############# Requirements ################## find_package(GSL REQUIRED) MESSAGE(STATUS "GSL include dirs: " ${GSL_INCLUDE_DIRS}) @@ -106,8 +109,6 @@ FILE(COPY ${HEADERS} DESTINATION ${PROJECT_BINARY_DIR}/include/catima) # the compiler used for C++ files MESSAGE( STATUS "CMAKE_CXX_COMPILER: " ${CMAKE_CXX_COMPILER} ) -# the compiler flags for compiling C++ sources -MESSAGE( STATUS "CMAKE_CXX_FLAGS: " ${CMAKE_CXX_FLAGS_RELEASE} ) ######## for python module if(PYTHON_MODULE) diff --git a/integrator.cpp b/integrator.cpp index 7b8809f..b0b7dcb 100644 --- a/integrator.cpp +++ b/integrator.cpp @@ -45,103 +45,5 @@ namespace catima{ return result; }; -//////////////////// GaussLegendre weights ///////////////////// -// order = 4 -template<> -std::array GaussLegendreIntegration<4>:: -x = {0.3399810435848562648026658,0.8611363115940525752239465}; - -template<> -std::array GaussLegendreIntegration<4>:: -w = {0.6521451548625461426269361,0.3478548451374538573730639}; - -// order = 6 -template<> -std::array GaussLegendreIntegration<6>:: -x = {0.2386191860831969086305017,0.6612093864662645136613996,0.9324695142031520278123016}; - -template<> -std::array GaussLegendreIntegration<6>:: -w = {0.4679139345726910473898703,0.3607615730481386075698335,0.1713244923791703450402961}; - -// order = 8 -template<> -std::array GaussLegendreIntegration<8>:: -x ={0.1834346424956498049394761,0.5255324099163289858177390,0.7966664774136267395915539,0.9602898564975362316835609}; - -template<> -std::array GaussLegendreIntegration<8>:: -w = {0.3626837833783619829651504,0.3137066458778872873379622,0.2223810344533744705443560,0.1012285362903762591525314}; - -// order = 10 -template<> -std::array GaussLegendreIntegration<10>:: - x = {0.1488743389816312108848260,0.4333953941292471907992659,0.6794095682990244062343274,0.8650633666889845107320967,0.9739065285171717200779640}; - -template<> -std::array GaussLegendreIntegration<10>:: - w = {0.2955242247147528701738930,0.2692667193099963550912269,0.2190863625159820439955349,0.1494513491505805931457763,0.0666713443086881375935688}; - -// order = 12 -template<> -std::array GaussLegendreIntegration<12>:: -x = {0.1252334085114689154724414,0.3678314989981801937526915,0.5873179542866174472967024,0.7699026741943046870368938,0.9041172563704748566784659,0.9815606342467192506905491}; - -template<> -std::array GaussLegendreIntegration<12>:: -w = {0.2491470458134027850005624,0.2334925365383548087608499,0.2031674267230659217490645,0.1600783285433462263346525,0.1069393259953184309602547,0.0471753363865118271946160}; - -// order = 14 -template<> -std::array GaussLegendreIntegration<14>:: -x = {0.1080549487073436620662447,0.3191123689278897604356718,0.5152486363581540919652907,0.6872929048116854701480198,0.8272013150697649931897947,0.9284348836635735173363911,0.9862838086968123388415973}; -template<> -std::array GaussLegendreIntegration<14>:: -w = {0.2152638534631577901958764,0.2051984637212956039659241,0.1855383974779378137417166,0.1572031671581935345696019,0.1215185706879031846894148,0.0801580871597602098056333,0.0351194603317518630318329}; - - -// order = 16 -template<> -std::array GaussLegendreIntegration<16>:: -x = {0.0950125098376374401853193,0.2816035507792589132304605,0.4580167776572273863424194,0.6178762444026437484466718,0.7554044083550030338951012,0.8656312023878317438804679,0.9445750230732325760779884,0.9894009349916499325961542}; - -template<> -std::array GaussLegendreIntegration<16>:: -w = {0.1894506104550684962853967,0.1826034150449235888667637,0.1691565193950025381893121,0.1495959888165767320815017,0.1246289712555338720524763,0.0951585116824927848099251,0.0622535239386478928628438,0.0271524594117540948517806}; - -// order = 20 -template<> -std::array GaussLegendreIntegration<20>:: -x = {0.0765265211334973337546404,0.2277858511416450780804962,0.3737060887154195606725482,0.5108670019508270980043641,0.6360536807265150254528367,0.7463319064601507926143051,0.8391169718222188233945291,0.9122344282513259058677524,0.9639719272779137912676661,0.9931285991850949247861224}; -template<> -std::array GaussLegendreIntegration<20>:: -w = {0.1527533871307258506980843,0.1491729864726037467878287,0.1420961093183820513292983,0.1316886384491766268984945,0.1181945319615184173123774,0.1019301198172404350367501,0.0832767415767047487247581,0.0626720483341090635695065,0.0406014298003869413310400,0.0176140071391521183118620}; - -// order = 32 -template<> -std::array GaussLegendreIntegration<32>:: -x = {0.0483076656877383162348126,0.1444719615827964934851864,0.2392873622521370745446032,0.3318686022821276497799168,0.4213512761306353453641194,0.5068999089322293900237475,0.5877157572407623290407455,0.6630442669302152009751152,0.7321821187402896803874267,0.7944837959679424069630973,0.8493676137325699701336930,0.8963211557660521239653072,0.9349060759377396891709191,0.9647622555875064307738119,0.9856115115452683354001750,0.9972638618494815635449811}; - -template<> -std::array GaussLegendreIntegration<32>:: -w = {0.0965400885147278005667648,0.0956387200792748594190820,0.0938443990808045656391802,0.0911738786957638847128686,0.0876520930044038111427715,0.0833119242269467552221991,0.0781938957870703064717409,0.0723457941088485062253994,0.0658222227763618468376501,0.0586840934785355471452836,0.0509980592623761761961632,0.0428358980222266806568786,0.0342738629130214331026877,0.0253920653092620594557526,0.0162743947309056706051706,0.0070186100094700966004071}; - -//order = 64 -template<> -std::array GaussLegendreIntegration<64>:: -x = {0.0243502926634244325089558,0.0729931217877990394495429,0.1214628192961205544703765,0.1696444204239928180373136,0.2174236437400070841496487,0.2646871622087674163739642,0.3113228719902109561575127,0.3572201583376681159504426,0.4022701579639916036957668,0.4463660172534640879849477,0.4894031457070529574785263,0.5312794640198945456580139,0.5718956462026340342838781,0.6111553551723932502488530,0.6489654712546573398577612,0.6852363130542332425635584,0.7198818501716108268489402,0.7528199072605318966118638,0.7839723589433414076102205,0.8132653151227975597419233,0.8406292962525803627516915,0.8659993981540928197607834,0.8893154459951141058534040,0.9105221370785028057563807,0.9295691721319395758214902,0.9464113748584028160624815,0.9610087996520537189186141,0.9733268277899109637418535,0.9833362538846259569312993,0.9910133714767443207393824,0.9963401167719552793469245,0.9993050417357721394569056}; - -template<> -std::array GaussLegendreIntegration<64>:: -w = {0.0486909570091397203833654,0.0485754674415034269347991,0.0483447622348029571697695,0.0479993885964583077281262,0.0475401657148303086622822,0.0469681828162100173253263,0.0462847965813144172959532,0.0454916279274181444797710,0.0445905581637565630601347,0.0435837245293234533768279,0.0424735151236535890073398,0.0412625632426235286101563,0.0399537411327203413866569,0.0385501531786156291289625,0.0370551285402400460404151,0.0354722132568823838106931,0.0338051618371416093915655,0.0320579283548515535854675,0.0302346570724024788679741,0.0283396726142594832275113,0.0263774697150546586716918,0.0243527025687108733381776,0.0222701738083832541592983,0.0201348231535302093723403,0.0179517157756973430850453,0.0157260304760247193219660,0.0134630478967186425980608,0.0111681394601311288185905,0.0088467598263639477230309,0.0065044579689783628561174,0.0041470332605624676352875,0.0017832807216964329472961}; - -//order = 128 -template<> -std::array GaussLegendreIntegration<128>:: -x = {0.0122236989606157641980521,0.0366637909687334933302153,0.0610819696041395681037870,0.0854636405045154986364980,0.1097942311276437466729747,0.1340591994611877851175753,0.1582440427142249339974755,0.1823343059853371824103826,0.2063155909020792171540580,0.2301735642266599864109866,0.2538939664226943208556180,0.2774626201779044028062316,0.3008654388776772026671541,0.3240884350244133751832523,0.3471177285976355084261628,0.3699395553498590266165917,0.3925402750332674427356482,0.4149063795522750154922739,0.4370245010371041629370429,0.4588814198335521954490891,0.4804640724041720258582757,0.5017595591361444642896063,0.5227551520511754784539479,0.5434383024128103634441936,0.5637966482266180839144308,0.5838180216287630895500389,0.6034904561585486242035732,0.6228021939105849107615396,0.6417416925623075571535249,0.6602976322726460521059468,0.6784589224477192593677557,0.6962147083695143323850866,0.7135543776835874133438599,0.7304675667419088064717369,0.7469441667970619811698824,0.7629743300440947227797691,0.7785484755064119668504941,0.7936572947621932902433329,0.8082917575079136601196422,0.8224431169556438424645942,0.8361029150609068471168753,0.8492629875779689691636001,0.8619154689395484605906323,0.8740527969580317986954180,0.8856677173453972174082924,0.8967532880491581843864474,0.9073028834017568139214859,0.9173101980809605370364836,0.9267692508789478433346245,0.9356743882779163757831268,0.9440202878302201821211114,0.9518019613412643862177963,0.9590147578536999280989185,0.9656543664319652686458290,0.9717168187471365809043384,0.9771984914639073871653744,0.9820961084357185360247656,0.9864067427245862088712355,0.9901278184917343833379303,0.9932571129002129353034372,0.9957927585349811868641612,0.9977332486255140198821574,0.9990774599773758950119878,0.9998248879471319144736081}; - -template<> -std::array GaussLegendreIntegration<128>:: -w = {0.0244461801962625182113259,0.0244315690978500450548486,0.0244023556338495820932980,0.0243585572646906258532685,0.0243002001679718653234426,0.0242273192228152481200933,0.0241399579890192849977167,0.0240381686810240526375873,0.0239220121367034556724504,0.0237915577810034006387807,0.0236468835844476151436514,0.0234880760165359131530253,0.0233152299940627601224157,0.0231284488243870278792979,0.0229278441436868469204110,0.0227135358502364613097126,0.0224856520327449668718246,0.0222443288937997651046291,0.0219897106684604914341221,0.0217219495380520753752610,0.0214412055392084601371119,0.0211476464682213485370195,0.0208414477807511491135839,0.0205227924869600694322850,0.0201918710421300411806732,0.0198488812328308622199444,0.0194940280587066028230219,0.0191275236099509454865185,0.0187495869405447086509195,0.0183604439373313432212893,0.0179603271850086859401969,0.0175494758271177046487069,0.0171281354231113768306810,0.0166965578015892045890915,0.0162550009097851870516575,0.0158037286593993468589656,0.0153430107688651440859909,0.0148731226021473142523855,0.0143943450041668461768239,0.0139069641329519852442880,0.0134112712886163323144890,0.0129075627392673472204428,0.0123961395439509229688217,0.0118773073727402795758911,0.0113513763240804166932817,0.0108186607395030762476596,0.0102794790158321571332153,0.0097341534150068058635483,0.0091830098716608743344787,0.0086263777986167497049788,0.0080645898904860579729286,0.0074979819256347286876720,0.0069268925668988135634267,0.0063516631617071887872143,0.0057726375428656985893346,0.0051901618326763302050708,0.0046045842567029551182905,0.0040162549837386423131943,0.0034255260409102157743378,0.0028327514714579910952857,0.0022382884309626187436221,0.0016425030186690295387909,0.0010458126793403487793129,0.0004493809602920903763943}; } diff --git a/integrator.h b/integrator.h index c396985..3bad088 100644 --- a/integrator.h +++ b/integrator.h @@ -46,7 +46,10 @@ class IntegratorGSL{ #endif }; -// built in integrator +template +struct GL_data{ +}; + template class GaussLegendreIntegration{ public: @@ -54,14 +57,10 @@ public: double integrate(F f, double a, double b) const; template double operator()(F f, double a, double b) const {return integrate(f, a, b);} - double get_w(int i) const {return w[i];} - double get_x(int i) const {return x[i];} + double w(int i) const {return GL_data::w()[i];} + double x(int i) const {return GL_data::x()[i];} int n() const {return order;} std::array get_points(double a = -1.0, double b = 1.0)const; - -public: - static std::array w; - static std::array x; }; template @@ -71,8 +70,9 @@ double GaussLegendreIntegration::integrate(F f, double a, double b) const double p = 0.5*(b-a); double q = 0.5*(b+a); - for(int i=0;i GaussLegendreIntegration::get_points(double a, int num = (order/2); for(int i=0;i< num;i++){ - points[num-i-1] = -p*x[i] + q; - points[num+i] = p*x[i] + q; + points[num-i-1] = -p*x(i) + q; + points[num+i] = p*x(i) + q; } return points; } +// order = 8 +template<> +struct GL_data<8>{ + static const std::array& x(){ + static const std::array _x = + {0.1834346424956498049394761,0.5255324099163289858177390,0.7966664774136267395915539,0.9602898564975362316835609}; + return _x; + } + static const std::array& w(){ + static const std::array _w = + {0.3626837833783619829651504,0.3137066458778872873379622,0.2223810344533744705443560,0.1012285362903762591525314}; + return _w; + } +}; + +// order = 10 +template<> +struct GL_data<10>{ + static const std::array& x(){ + static const std::array _x = + {0.1488743389816312108848260,0.4333953941292471907992659,0.6794095682990244062343274,0.8650633666889845107320967,0.9739065285171717200779640}; + return _x; + } + static const std::array& w(){ + static const std::array _w = + {0.2955242247147528701738930,0.2692667193099963550912269,0.2190863625159820439955349,0.1494513491505805931457763,0.0666713443086881375935688}; + return _w; + } +}; + +// order = 12 +template<> +struct GL_data<12>{ + static const std::array& x(){ + static const std::array _x = + {0.1252334085114689154724414,0.3678314989981801937526915,0.5873179542866174472967024,0.7699026741943046870368938,0.9041172563704748566784659,0.9815606342467192506905491}; + return _x; + } + + static const std::array& w(){ + static const std::array _w = + {0.2491470458134027850005624,0.2334925365383548087608499,0.2031674267230659217490645,0.1600783285433462263346525,0.1069393259953184309602547,0.0471753363865118271946160}; + return _w; + } +}; + +// order = 14 +template<> +struct GL_data<14>{ + static const std::array& x(){ + static const std::array _x = + {0.1080549487073436620662447,0.3191123689278897604356718,0.5152486363581540919652907,0.6872929048116854701480198,0.8272013150697649931897947,0.9284348836635735173363911,0.9862838086968123388415973}; + return _x; + } + + static const std::array& w(){ + static const std::array _w = + {0.2152638534631577901958764,0.2051984637212956039659241,0.1855383974779378137417166,0.1572031671581935345696019,0.1215185706879031846894148,0.0801580871597602098056333,0.0351194603317518630318329}; + return _w; + } +}; + +// order = 16 +template<> +struct GL_data<16>{ + static const std::array& x(){ + static const std::array _x = + {0.0950125098376374401853193,0.2816035507792589132304605,0.4580167776572273863424194,0.6178762444026437484466718,0.7554044083550030338951012,0.8656312023878317438804679,0.9445750230732325760779884,0.9894009349916499325961542}; + return _x; + } + + static const std::array& w(){ + static const std::array _w = + {0.1894506104550684962853967,0.1826034150449235888667637,0.1691565193950025381893121,0.1495959888165767320815017,0.1246289712555338720524763,0.0951585116824927848099251,0.0622535239386478928628438,0.0271524594117540948517806}; + return _w; + } +}; + #ifdef GSL_INTEGRATION using integrator_type = IntegratorGSL; #else diff --git a/tests/test.py b/tests/test.py index 7978ad6..da604c6 100644 --- a/tests/test.py +++ b/tests/test.py @@ -169,6 +169,23 @@ class TestStructures(unittest.TestCase): res = catima.calculate(p,water) res2 = catima.dedx_from_range(p,water) self.assertAlmostEqual(res.dEdxi,res2,3) + + def test_config(self): + water = catima.get_material(catima.material.WATER) + water.density(1.0) + water.thickness(1.0) + p = catima.Projectile(1,1) + conf = catima.Config() + conf.dedx_straggling(catima.omega_type.bohr) + conf2 = catima.Config() + conf2.dedx_straggling(catima.omega_type.atima) + p(1000) + res = catima.calculate(p,water,config=conf) + res2 = catima.calculate(p,water,config=conf2) + self.assertAlmostEqual(res.dEdxi,res2.dEdxi,delta=1e-6) + self.assertNotAlmostEqual(res.sigma_E,res2.sigma_E,delta=1e-4) + self.assertNotAlmostEqual(res.sigma_r,res2.sigma_r,delta=1e-4) + def test_eout(self): graphite = catima.get_material(6) diff --git a/tests/test_dedx_range.cpp b/tests/test_dedx_range.cpp index d6d6df4..6c657be 100644 --- a/tests/test_dedx_range.cpp +++ b/tests/test_dedx_range.cpp @@ -9,8 +9,8 @@ using std::cout; using std::endl; -double logEmax = catima::logEmax; -double logEmin = catima::logEmin; +double logEmax = catima::logEmax-0.01; +double logEmin = catima::logEmin+0.001; template struct EnergyTable{ @@ -106,7 +106,7 @@ int main(){ std::vector> projectiles = {{1,1},{4,2},{12,6},{238,92}}; bool fwrite = false; - double eps = 0.01; + double eps = 0.005; comp_dedx(p,c_target,eps,fwrite,"dedx_p_c.dat"); comp_dedx(p,h_target,eps,fwrite,"dedx_p_h.dat"); comp_dedx(p,water,eps,fwrite,"dedx_p_water.dat");