From f9f3aa6e96e9c7827bef1f449fa2ce3d86505a23 Mon Sep 17 00:00:00 2001 From: Luca Fedeli Date: Mon, 21 Sep 2020 12:42:02 +0200 Subject: Coupling WarpX with an ✨improved✨ version of the QED library (#1198) MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit * Initial work to couple improved QED module to WarpX * WIP to couple with WarpX the new QED library * Continuing work to couple the new version of the QED library with WarpX * progress towards completing coupling with new version of QED library * WarpX coupled with new version of QED library * default behavior is to display table generation progress * some host device functions are now device only * fixed bug * bugfixing * updating tests * updated test * updated test * added initial version of tests (not working) * added check and updated a comment * fixed bug * added inputfiles and analysis script for new BW tests * test for BW process are ready * modified test * make lgtm happy * removed TABs * initial work to add QS tests (not working) * removed old tests * fixed bug in script * changed position of evolution of optical depth * progress with QSR tests * improved test * very low energy photons are always eliminated * added tests to regression suite * improved test * improved tests * removed redundant parameter * removed trailing white space * updated documentation * fix lgtm warnings * fixed missing check on chi parameter * fixed missing check on chi parameter & bugfixing * improved comments * increased tolerance in tests * updated units in test * now test succeds if the error is extremely small * updated checksums * fixed bug * fixed some unused or uninitialized variables warnings * now using ignore_unused instead of commenting out some variables * fixed warnings * partial fix of a test * fixed test * fixed test * added checksums * fixed tests * fixed benchmark for qed_schwinger2 * removed checksums for tests which do no exist anymore * fixed checksums for several qed tests * fixed checksums for several qed tests * fixed checksums * removed unwanted checksum * fixed checksum * removed files which should have been deleted * add some const * [skip ci] added some docstrings and some const * Update Source/Particles/ElementaryProcess/QEDInternals/BreitWheelerEngineWrapper.H Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/ElementaryProcess/QEDInternals/BreitWheelerEngineWrapper.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/ElementaryProcess/QEDInternals/QuantumSyncEngineWrapper.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * added some docstrings and some const * replaced ManagedVectors with DeviceVectors * Update Source/Particles/ElementaryProcess/QEDInternals/QedWrapperCommons.H Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * added some const * removed unwanted assert * updated comment * changed position of GPU synchronization directive * Update Docs/source/running_cpp/parameters.rst Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Examples/Modules/qed/quantum_synchrotron/analysis.py Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Examples/Modules/qed/quantum_synchrotron/analysis.py Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Examples/Modules/qed/breit_wheeler/analysis.py Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Examples/Modules/qed/breit_wheeler/analysis.py Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * add do_plot option to some analysis scripts * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * uncomment a line * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * simplified input scripts for BW tests * simplified input scripts for QS tests * removed unwanted files * simplified analysis script * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * reverted modification to schwinger analysis script * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * remove outdated comment * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Source/Particles/MultiParticleContainer.cpp Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * fix warnings * made test more robust * reset benchmark for qed_breit_wheeler_2d * fixed bug in test * make test more robust * made test more robust * Update Examples/Modules/qed/quantum_synchrotron/analysis.py Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update Examples/Modules/qed/quantum_synchrotron/analysis.py Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> * Update run_test.sh Co-authored-by: Axel Huebl Co-authored-by: Tools Co-authored-by: NeilZaim <49716072+NeilZaim@users.noreply.github.com> Co-authored-by: Axel Huebl --- Source/Particles/MultiParticleContainer.cpp | 185 ++++++++++++++-------------- 1 file changed, 95 insertions(+), 90 deletions(-) (limited to 'Source/Particles/MultiParticleContainer.cpp') diff --git a/Source/Particles/MultiParticleContainer.cpp b/Source/Particles/MultiParticleContainer.cpp index a866ae3a0..cb47755f7 100644 --- a/Source/Particles/MultiParticleContainer.cpp +++ b/Source/Particles/MultiParticleContainer.cpp @@ -24,7 +24,7 @@ #include #include #include - +#include using namespace amrex; @@ -756,7 +756,9 @@ void MultiParticleContainer::InitQuantumSync () //If specified, use a user-defined energy threshold for photon creaction ParticleReal temp; + constexpr auto mec2 = PhysConst::c * PhysConst::c * PhysConst::m_e; if(pp.query("photon_creation_energy_threshold", temp)){ + temp *= mec2; m_quantum_sync_photon_creation_energy_threshold = temp; } else{ @@ -764,6 +766,13 @@ void MultiParticleContainer::InitQuantumSync () " for photon energy creaction threshold \n" ; } + // qs_minimum_chi_part is the minimum chi parameter to be + // considered for Synchrotron emission. If a lepton has chi < chi_min, + // the optical depth is not evolved and photon generation is ignored + amrex::Real qs_minimum_chi_part; + pp.get("chi_min", qs_minimum_chi_part); + + pp.query("lookup_table_mode", lookup_table_mode); if(lookup_table_mode.empty()){ amrex::Abort("Quantum Synchrotron table mode should be provided"); @@ -787,11 +796,12 @@ void MultiParticleContainer::InitQuantumSync () Vector table_data; ParallelDescriptor::ReadAndBcastFile(load_table_name, table_data); ParallelDescriptor::Barrier(); - m_shr_p_qs_engine->init_lookup_tables_from_raw_data(table_data); + m_shr_p_qs_engine->init_lookup_tables_from_raw_data(table_data, + qs_minimum_chi_part); } - else if(lookup_table_mode == "dummy_builtin"){ - amrex::Print() << "Built-in Quantum Synchrotron dummy table will be used. \n" ; - m_shr_p_qs_engine->init_dummy_tables(); + else if(lookup_table_mode == "builtin"){ + amrex::Print() << "Built-in Quantum Synchrotron table will be used. \n" ; + m_shr_p_qs_engine->init_builtin_tables(qs_minimum_chi_part); } else{ amrex::Abort("Unknown Quantum Synchrotron table mode"); @@ -806,6 +816,14 @@ void MultiParticleContainer::InitBreitWheeler () { std::string lookup_table_mode; ParmParse pp("qed_bw"); + + // bw_minimum_chi_phot is the minimum chi parameter to be + // considered for pair production. If a photon has chi < chi_min, + // the optical depth is not evolved and photon generation is ignored + amrex::Real bw_minimum_chi_part; + if(!pp.query("chi_min", bw_minimum_chi_part)) + amrex::Abort("qed_bw.chi_min should be provided!"); + pp.query("lookup_table_mode", lookup_table_mode); if(lookup_table_mode.empty()){ amrex::Abort("Breit Wheeler table mode should be provided"); @@ -829,11 +847,12 @@ void MultiParticleContainer::InitBreitWheeler () Vector table_data; ParallelDescriptor::ReadAndBcastFile(load_table_name, table_data); ParallelDescriptor::Barrier(); - m_shr_p_bw_engine->init_lookup_tables_from_raw_data(table_data); + m_shr_p_bw_engine->init_lookup_tables_from_raw_data( + table_data, bw_minimum_chi_part); } - else if(lookup_table_mode == "dummy_builtin"){ - amrex::Print() << "Built-in Breit Wheeler dummy table will be used. \n" ; - m_shr_p_bw_engine->init_dummy_tables(); + else if(lookup_table_mode == "builtin"){ + amrex::Print() << "Built-in Breit Wheeler table will be used. \n" ; + m_shr_p_bw_engine->init_builtin_tables(bw_minimum_chi_part); } else{ amrex::Abort("Unknown Breit Wheeler table mode"); @@ -853,15 +872,14 @@ MultiParticleContainer::QuantumSyncGenerateTable () if(table_name.empty()) amrex::Abort("qed_qs.save_table_in should be provided!"); - if(ParallelDescriptor::IOProcessor()){ - PicsarQuantumSynchrotronCtrl ctrl; - int t_int; + // qs_minimum_chi_part is the minimum chi parameter to be + // considered for Synchrotron emission. If a lepton has chi < chi_min, + // the optical depth is not evolved and photon generation is ignored + amrex::Real qs_minimum_chi_part; + pp.get("chi_min", qs_minimum_chi_part); - // Engine paramenter: chi_part_min is the minium chi parameter to be - // considered by the engine. If a lepton has chi < chi_part_min, - // the optical depth is not evolved and photon generation is ignored - if(!pp.query("chi_min", ctrl.chi_part_min)) - amrex::Abort("qed_qs.chi_min should be provided!"); + if(ParallelDescriptor::IOProcessor()){ + PicsarQuantumSyncCtrl ctrl; //==Table parameters== @@ -869,20 +887,16 @@ MultiParticleContainer::QuantumSyncGenerateTable () //These parameters are used to pre-compute a function //which appears in the evolution of the optical depth - //Minimun chi for the table. If a lepton has chi < chi_part_tdndt_min, - //chi is considered as it were equal to chi_part_tdndt_min - if(!pp.query("tab_dndt_chi_min", ctrl.chi_part_tdndt_min)) - amrex::Abort("qed_qs.tab_dndt_chi_min should be provided!"); + //Minimun chi for the table. If a lepton has chi < tab_dndt_chi_min, + //chi is considered as if it were equal to tab_dndt_chi_min + pp.get("tab_dndt_chi_min", ctrl.dndt_params.chi_part_min); - //Maximum chi for the table. If a lepton has chi > chi_part_tdndt_max, - //chi is considered as it were equal to chi_part_tdndt_max - if(!pp.query("tab_dndt_chi_max", ctrl.chi_part_tdndt_max)) - amrex::Abort("qed_qs.tab_dndt_chi_max should be provided!"); + //Maximum chi for the table. If a lepton has chi > tab_dndt_chi_max, + //chi is considered as if it were equal to tab_dndt_chi_max + pp.get("tab_dndt_chi_max", ctrl.dndt_params.chi_part_max); //How many points should be used for chi in the table - if(!pp.query("tab_dndt_how_many", t_int)) - amrex::Abort("qed_qs.tab_dndt_how_many should be provided!"); - ctrl.chi_part_tdndt_how_many = t_int; + pp.get("tab_dndt_how_many", ctrl.dndt_params.chi_part_how_many); //------ //--- sub-table 2 (2D) @@ -890,32 +904,31 @@ MultiParticleContainer::QuantumSyncGenerateTable () //which is used to extract the properties of the generated //photons. - //Minimun chi for the table. If a lepton has chi < chi_part_tem_min, - //chi is considered as it were equal to chi_part_tem_min - if(!pp.query("tab_em_chi_min", ctrl.chi_part_tem_min)) - amrex::Abort("qed_qs.tab_em_chi_min should be provided!"); + //Minimun chi for the table. If a lepton has chi < tab_em_chi_min, + //chi is considered as if it were equal to tab_em_chi_min + pp.get("tab_em_chi_min", ctrl.phot_em_params.chi_part_min); - //Maximum chi for the table. If a lepton has chi > chi_part_tem_max, - //chi is considered as it were equal to chi_part_tem_max - if(!pp.query("tab_em_chi_max", ctrl.chi_part_tem_max)) - amrex::Abort("qed_qs.tab_em_chi_max should be provided!"); + //Maximum chi for the table. If a lepton has chi > tab_em_chi_max, + //chi is considered as if it were equal to tab_em_chi_max + pp.get("tab_em_chi_max", ctrl.phot_em_params.chi_part_max); //How many points should be used for chi in the table - if(!pp.query("tab_em_chi_how_many", t_int)) - amrex::Abort("qed_qs.tab_em_chi_how_many should be provided!"); - ctrl.chi_part_tem_how_many = t_int; - - //The other axis of the table is a cumulative probability distribution - //(corresponding to different energies of the generated particles) - //This parameter is the number of different points to consider - if(!pp.query("tab_em_prob_how_many", t_int)) - amrex::Abort("qed_qs.tab_em_prob_how_many should be provided!"); - ctrl.prob_tem_how_many = t_int; + pp.get("tab_em_chi_how_many", ctrl.phot_em_params.chi_part_how_many); + + //The other axis of the table is the ratio between the quantum + //parameter of the emitted photon and the quantum parameter of the + //lepton. This parameter is the minimum ratio to consider for the table. + pp.get("tab_em_frac_min", ctrl.phot_em_params.frac_min); + + //This parameter is the number of different points to consider for the second + //axis + pp.get("tab_em_frac_how_many", ctrl.phot_em_params.frac_how_many); //==================== - m_shr_p_qs_engine->compute_lookup_tables(ctrl); + m_shr_p_qs_engine->compute_lookup_tables(ctrl, qs_minimum_chi_part); + const auto data = m_shr_p_qs_engine->export_lookup_tables_data(); WarpXUtilIO::WriteBinaryDataOnFile(table_name, - m_shr_p_qs_engine->export_lookup_tables_data()); + Vector{data.begin(), data.end()}); } ParallelDescriptor::Barrier(); @@ -926,7 +939,8 @@ MultiParticleContainer::QuantumSyncGenerateTable () //No need to initialize from raw data for the processor that //has just generated the table if(!ParallelDescriptor::IOProcessor()){ - m_shr_p_qs_engine->init_lookup_tables_from_raw_data(table_data); + m_shr_p_qs_engine->init_lookup_tables_from_raw_data( + table_data, qs_minimum_chi_part); } } @@ -939,15 +953,14 @@ MultiParticleContainer::BreitWheelerGenerateTable () if(table_name.empty()) amrex::Abort("qed_bw.save_table_in should be provided!"); + // bw_minimum_chi_phot is the minimum chi parameter to be + // considered for pair production. If a photon has chi < chi_min, + // the optical depth is not evolved and photon generation is ignored + amrex::Real bw_minimum_chi_part; + pp.get("chi_min", bw_minimum_chi_part); + if(ParallelDescriptor::IOProcessor()){ PicsarBreitWheelerCtrl ctrl; - int t_int; - - // Engine paramenter: chi_phot_min is the minium chi parameter to be - // considered by the engine. If a photon has chi < chi_phot_min, - // the optical depth is not evolved and pair generation is ignored - if(!pp.query("chi_min", ctrl.chi_phot_min)) - amrex::Abort("qed_bw.chi_min should be provided!"); //==Table parameters== @@ -955,20 +968,16 @@ MultiParticleContainer::BreitWheelerGenerateTable () //These parameters are used to pre-compute a function //which appears in the evolution of the optical depth - //Minimun chi for the table. If a photon has chi < chi_phot_tdndt_min, + //Minimun chi for the table. If a photon has chi < tab_dndt_chi_min, //an analytical approximation is used. - if(!pp.query("tab_dndt_chi_min", ctrl.chi_phot_tdndt_min)) - amrex::Abort("qed_bw.tab_dndt_chi_min should be provided!"); + pp.get("tab_dndt_chi_min", ctrl.dndt_params.chi_phot_min); - //Maximum chi for the table. If a photon has chi > chi_phot_tdndt_min, + //Maximum chi for the table. If a photon has chi > tab_dndt_chi_max, //an analytical approximation is used. - if(!pp.query("tab_dndt_chi_max", ctrl.chi_phot_tdndt_max)) - amrex::Abort("qed_bw.tab_dndt_chi_max should be provided!"); + pp.get("tab_dndt_chi_max", ctrl.dndt_params.chi_phot_max); //How many points should be used for chi in the table - if(!pp.query("tab_dndt_how_many", t_int)) - amrex::Abort("qed_bw.tab_dndt_how_many should be provided!"); - ctrl.chi_phot_tdndt_how_many = t_int; + pp.get("tab_dndt_how_many", ctrl.dndt_params.chi_phot_how_many); //------ //--- sub-table 2 (2D) @@ -976,32 +985,27 @@ MultiParticleContainer::BreitWheelerGenerateTable () //which is used to extract the properties of the generated //particles. - //Minimun chi for the table. If a photon has chi < chi_phot_tpair_min + //Minimun chi for the table. If a photon has chi < tab_pair_chi_min //chi is considered as it were equal to chi_phot_tpair_min - if(!pp.query("tab_pair_chi_min", ctrl.chi_phot_tpair_min)) - amrex::Abort("qed_bw.tab_pair_chi_min should be provided!"); + pp.get("tab_pair_chi_min", ctrl.pair_prod_params.chi_phot_min); - //Maximum chi for the table. If a photon has chi > chi_phot_tpair_max + //Maximum chi for the table. If a photon has chi > tab_pair_chi_max //chi is considered as it were equal to chi_phot_tpair_max - if(!pp.query("tab_pair_chi_max", ctrl.chi_phot_tpair_max)) - amrex::Abort("qed_bw.tab_pair_chi_max should be provided!"); + pp.get("tab_pair_chi_max", ctrl.pair_prod_params.chi_phot_max); //How many points should be used for chi in the table - if(!pp.query("tab_pair_chi_how_many", t_int)) - amrex::Abort("qed_bw.tab_pair_chi_how_many should be provided!"); - ctrl.chi_phot_tpair_how_many = t_int; + pp.get("tab_pair_chi_how_many", ctrl.pair_prod_params.chi_phot_how_many); //The other axis of the table is the fraction of the initial energy //'taken away' by the most energetic particle of the pair. //This parameter is the number of different fractions to consider - if(!pp.query("tab_pair_frac_how_many", t_int)) - amrex::Abort("qed_bw.tab_pair_frac_how_many should be provided!"); - ctrl.chi_frac_tpair_how_many = t_int; + pp.get("tab_pair_frac_how_many", ctrl.pair_prod_params.frac_how_many); //==================== - m_shr_p_bw_engine->compute_lookup_tables(ctrl); + m_shr_p_bw_engine->compute_lookup_tables(ctrl, bw_minimum_chi_part); + const auto data = m_shr_p_bw_engine->export_lookup_tables_data(); WarpXUtilIO::WriteBinaryDataOnFile(table_name, - m_shr_p_bw_engine->export_lookup_tables_data()); + Vector{data.begin(), data.end()}); } ParallelDescriptor::Barrier(); @@ -1012,7 +1016,8 @@ MultiParticleContainer::BreitWheelerGenerateTable () //No need to initialize from raw data for the processor that //has just generated the table if(!ParallelDescriptor::IOProcessor()){ - m_shr_p_bw_engine->init_lookup_tables_from_raw_data(table_data); + m_shr_p_bw_engine->init_lookup_tables_from_raw_data( + table_data, bw_minimum_chi_part); } } @@ -1039,21 +1044,20 @@ MultiParticleContainer::doQEDSchwinger () amrex::Abort("Schwinger process not implemented in rz geometry"); #endif -#ifdef AMREX_USE_FLOAT - amrex::Abort("Schwinger process not implemented in single precision"); -#endif - -// Get cell volume multiplied by temporal step. In 2D the transverse size is +// Get cell volume. In 2D the transverse size is // chosen by the user in the input file. amrex::Geometry const & geom = warpx.Geom(level_0); #if (AMREX_SPACEDIM == 2) - const auto dVdt = geom.CellSize(0) * geom.CellSize(1) - * m_qed_schwinger_y_size * warpx.getdt(level_0); + const auto dV = geom.CellSize(0) * geom.CellSize(1) + * m_qed_schwinger_y_size; #elif (AMREX_SPACEDIM == 3) - const auto dVdt = geom.CellSize(0) * geom.CellSize(1) - * geom.CellSize(2) * warpx.getdt(level_0); + const auto dV = geom.CellSize(0) * geom.CellSize(1) + * geom.CellSize(2); #endif + // Get the temporal step + const auto dt = warpx.getdt(level_0); + auto& pc_product_ele = allcontainers[m_qed_schwinger_ele_product]; auto& pc_product_pos = @@ -1100,7 +1104,8 @@ MultiParticleContainer::doQEDSchwinger () const auto np_pos_dst = dst_pos_tile.numParticles(); const auto Filter = SchwingerFilterFunc{ - m_qed_schwinger_threshold_poisson_gaussian,dVdt}; + m_qed_schwinger_threshold_poisson_gaussian, + dV, dt}; const SmartCreateFactory create_factory_ele(*pc_product_ele); const SmartCreateFactory create_factory_pos(*pc_product_pos); -- cgit v1.2.3