diff options
Diffstat (limited to 'Source/Particles/MultiParticleContainer.cpp')
-rw-r--r-- | Source/Particles/MultiParticleContainer.cpp | 185 |
1 files changed, 95 insertions, 90 deletions
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 <limits> #include <algorithm> #include <string> - +#include <vector> 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<char> 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<char> 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<char>{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<char>{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); |