From a08c76c7857d619f39d51053ed97129172ff8799 Mon Sep 17 00:00:00 2001 From: NeilZaim <49716072+NeilZaim@users.noreply.github.com> Date: Mon, 11 May 2020 17:53:07 +0200 Subject: Add Schwinger process (#784) * Initial work to add back QED particle generation * Work in progress: port old QED routines * Add two distinct CopyFuncs * modified getMFItInfo and CopyFunc (not working) * bugfixing & work to add back QED particle creation routines * bugfixing * added back quantum photon emission * bugfixing * bugfixing * added back pair generation (still some bugs in photon emission) * removed unwanted check * bugfixing * bugfixing * bugfixing * Moved QED folder * added comments + some refactoring * added comments * remove some virtual functions to make lgtm happy * updated tests * added PhysicalParticleType * bugfixing * added copyright * improved comments * improved comments * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * moved inclusion of QEDInternals folder * moved some inclusion directives between Make files * moved some inclusion directives between Make files (forgot to add a file) * Update Source/Particles/PhysicalParticleContainer.H Co-Authored-By: Axel Huebl * Update Source/Particles/PhysicalParticleContainer.H Co-Authored-By: Axel Huebl * Update Source/Particles/PhysicalParticleContainer.H Co-Authored-By: Axel Huebl * Update Source/Particles/PhysicalParticleContainer.H Co-Authored-By: Axel Huebl * Update Source/Particles/PhysicalParticleContainer.H Co-Authored-By: Axel Huebl * Update Source/Particles/MultiParticleContainer.cpp Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: Axel Huebl * Update Source/Particles/MultiParticleContainer.cpp Co-Authored-By: MaxThevenet * Update Source/Particles/MultiParticleContainer.cpp Co-Authored-By: MaxThevenet * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: Axel Huebl * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: Axel Huebl * Update Source/Particles/MultiParticleContainer.cpp Co-Authored-By: MaxThevenet * Update Source/Particles/MultiParticleContainer.cpp Co-Authored-By: MaxThevenet * Update Source/Particles/MultiParticleContainer.cpp Co-Authored-By: MaxThevenet * Update Source/Particles/MultiParticleContainer.H Co-Authored-By: Axel Huebl * corrected alignment * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: Axel Huebl * removed some unnecessary amrex:: * add missing comment * Replaced BL_PROFILE with WARPX_PROFILE * bugfixing and making some variables const * removed some moves * removed some moves * started to change tau into optical_depth_BW or optical_depth_QSR * Using initialization policy to initialize optical depth * bugfixing * forgot to add a file * fixed bug * Revert "fixed bug" This reverts commit a3fb98d10cc30327635aeaa71451a05ca2229ff4. * Define doQEDSchwinger function * Read input parameters for Schwinger process * Update Source/Particles/ElementaryProcess/QEDPairGeneration.H Co-Authored-By: MaxThevenet * Update Source/Particles/MultiParticleContainer.cpp Co-Authored-By: MaxThevenet * Added doQEDEvents to OneStep_sub1 * add a bunch of const * add _rt suffix * Update Source/Particles/MultiParticleContainer.cpp Co-Authored-By: MaxThevenet * added path to included files * Introduced a templated AmIA function using physical_species * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: MaxThevenet * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: MaxThevenet * Update Source/Particles/ElementaryProcess/QEDPhotonEmission.H Co-Authored-By: MaxThevenet * added paths to included headers * updated documentation * updated examples * bugfixing * bugfixing * fixing examples * fixed example * fixed example * correct a misprint in error message * fixed issue related to 1./mass for photons * Add skeleton of doQEDSchwinger function * Linked to PICSAR to calculate pair production rate * Added missing header file * fix conflicts in QEDPhotonEmission.H * Written first version of FilterCreateTransformFromFAB.H * Cleanup some useless comments * Update Source/Particles/ElementaryProcess/QEDInternals/SchwingerProcessWrapper.H Co-Authored-By: Luca Fedeli * Update Source/Particles/ParticleCreation/FilterCreateTransformFromFAB.H Co-Authored-By: Luca Fedeli * Update Source/Particles/ParticleCreation/FilterCreateTransformFromFAB.H Co-Authored-By: Luca Fedeli * Minor revisions and improvements * Remove trailing white spaces * Write filter function * Remove print * Some debugging and cleaning * Write tranform function * remove transform_dummy * Added some tests * Remove EOL whitespaces * update prepare_file_travis.py * Should fix error in automated tests * Make the tests run in parallel * Put path relative to Source/ in includes * update include path * Actually resolve conflicts * put dVdt and weight_index as members of filter and transform functions * a bit of debugging on GPU * Add comments and documentation * Fix typos * Add assert for single precision * Update Docs/source/running_cpp/parameters.rst Co-Authored-By: Luca Fedeli * Update Source/Particles/MultiParticleContainer.H Co-Authored-By: Luca Fedeli * Update Source/Particles/ParticleCreation/FilterCreateTransformFromFAB.H Co-Authored-By: Luca Fedeli * Add profiler and minor modifications * Fix typo * Update asserts so that module works with momentum conserving algo * update assert if particles are created in analysis script * Apply suggestions from code review Co-Authored-By: MaxThevenet * Add comments to FilterCopyTransform and FilterCreateTransform functions * Add const and comment * Update Source/Particles/ParticleCreation/FilterCreateTransformFromFAB.H Co-authored-by: MaxThevenet * Update test for new diags and remove AllContainerType * Update Regression/WarpX-tests.ini * update test * update tests * Split Schwinger test into 4 separate tests * Change name of analysis script * update tests * Combine all analysis scripts into a single one Co-authored-by: Luca Fedeli Co-authored-by: Luca Fedeli Co-authored-by: Axel Huebl Co-authored-by: MaxThevenet Co-authored-by: Neil --- Source/Particles/MultiParticleContainer.cpp | 150 +++++++++++++++++++++++++++- 1 file changed, 148 insertions(+), 2 deletions(-) (limited to 'Source/Particles/MultiParticleContainer.cpp') diff --git a/Source/Particles/MultiParticleContainer.cpp b/Source/Particles/MultiParticleContainer.cpp index a3b5d256e..123771c84 100644 --- a/Source/Particles/MultiParticleContainer.cpp +++ b/Source/Particles/MultiParticleContainer.cpp @@ -1,7 +1,8 @@ /* Copyright 2019-2020 Andrew Myers, Ann Almgren, Axel Huebl * David Grote, Jean-Luc Vay, Luca Fedeli - * Mathieu Lobet, Maxence Thevenet, Remi Lehe - * Revathi Jambunathan, Weiqun Zhang, Yinjian Zhao + * Mathieu Lobet, Maxence Thevenet, Neil Zaim + * Remi Lehe, Revathi Jambunathan, Weiqun Zhang + * Yinjian Zhao * * * This file is part of WarpX. @@ -12,6 +13,11 @@ #include "SpeciesPhysicalProperties.H" #include "Utils/WarpXUtil.H" #include "WarpX.H" +#ifdef WARPX_QED + #include "Particles/ElementaryProcess/QEDInternals/SchwingerProcessWrapper.H" + #include "Particles/ElementaryProcess/QEDSchwingerProcess.H" + #include "Particles/ParticleCreation/FilterCreateTransformFromFAB.H" +#endif #include @@ -258,6 +264,21 @@ MultiParticleContainer::ReadParameters () AMREX_ALWAYS_ASSERT(lasers_names.size() == nlasers); } +#ifdef WARPX_QED + ParmParse ppw("warpx"); + ppw.query("do_qed_schwinger", m_do_qed_schwinger); + + if (m_do_qed_schwinger) { + ParmParse ppq("qed_schwinger"); + ppq.get("ele_product_species", m_qed_schwinger_ele_product_name); + ppq.get("pos_product_species", m_qed_schwinger_pos_product_name); +#if (AMREX_SPACEDIM == 2) + ppq.get("y_size",m_qed_schwinger_y_size); +#endif + ppq.query("threshold_poisson_gaussian", + m_qed_schwinger_threshold_poisson_gaussian); + } +#endif initialized = true; } } @@ -576,6 +597,14 @@ MultiParticleContainer::mapSpeciesProduct () } +#ifdef WARPX_QED + if (m_do_qed_schwinger) { + m_qed_schwinger_ele_product = + getSpeciesID(m_qed_schwinger_ele_product_name); + m_qed_schwinger_pos_product = + getSpeciesID(m_qed_schwinger_pos_product_name); + } +#endif } /* \brief Given a species name, return its ID. @@ -991,6 +1020,112 @@ MultiParticleContainer::BreitWheelerGenerateTable () } } +void +MultiParticleContainer::doQEDSchwinger () +{ + WARPX_PROFILE("MPC::doQEDSchwinger"); + + if (!m_do_qed_schwinger) {return;} + + auto & warpx = WarpX::GetInstance(); + + AMREX_ALWAYS_ASSERT_WITH_MESSAGE(warpx.do_nodal || + warpx.field_gathering_algo == GatheringAlgo::MomentumConserving, + "ERROR: Schwinger process only implemented for warpx.do_nodal = 1" + "or algo.field_gathering = momentum-conserving"); + + const int level_0 = 0; + + AMREX_ALWAYS_ASSERT_WITH_MESSAGE(warpx.maxLevel() == level_0, + "ERROR: Schwinger process not implemented with mesh refinement"); + +#ifdef WARPX_DIM_RZ + 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 +// chosen by the user in the input file. + amrex::Geometry const & geom = warpx.Geom(level_0); + auto domain_box = geom.Domain(); +#if (AMREX_SPACEDIM == 2) + const auto dVdt = geom.CellSize(0) * geom.CellSize(1) + * m_qed_schwinger_y_size * warpx.getdt(level_0); +#elif (AMREX_SPACEDIM == 3) + const auto dVdt = geom.CellSize(0) * geom.CellSize(1) + * geom.CellSize(2) * warpx.getdt(level_0); +#endif + + auto& pc_product_ele = + allcontainers[m_qed_schwinger_ele_product]; + auto& pc_product_pos = + allcontainers[m_qed_schwinger_pos_product]; + + const MultiFab & Ex = warpx.getEfield(level_0,0); + const MultiFab & Ey = warpx.getEfield(level_0,1); + const MultiFab & Ez = warpx.getEfield(level_0,2); + const MultiFab & Bx = warpx.getBfield(level_0,0); + const MultiFab & By = warpx.getBfield(level_0,1); + const MultiFab & Bz = warpx.getBfield(level_0,2); + + MFItInfo info; + if (TilingIfNotGPU()) { + info.EnableTiling(); + } +#ifdef _OPENMP + info.SetDynamic(WarpX::do_dynamic_scheduling); +#pragma omp parallel +#endif + + for (MFIter mfi(Ex, info); mfi.isValid(); ++mfi ) + { + // Make the box cell centered to avoid creating particles twice on the tile edges + const Box& box = enclosedCells(mfi.tilebox()); + + const auto& arrEx = Ex[mfi].array(); + const auto& arrEy = Ey[mfi].array(); + const auto& arrEz = Ez[mfi].array(); + const auto& arrBx = Bx[mfi].array(); + const auto& arrBy = By[mfi].array(); + const auto& arrBz = Bz[mfi].array(); + + const Array4 array_EMFAB [] = {arrEx,arrEy,arrEz, + arrBx,arrBy,arrBz}; + + pc_product_ele->defineAllParticleTiles(); + pc_product_pos->defineAllParticleTiles(); + + auto& dst_ele_tile = pc_product_ele->ParticlesAt(level_0, mfi); + auto& dst_pos_tile = pc_product_pos->ParticlesAt(level_0, mfi); + + const auto np_ele_dst = dst_ele_tile.numParticles(); + const auto np_pos_dst = dst_pos_tile.numParticles(); + + const auto Filter = SchwingerFilterFunc{ + m_qed_schwinger_threshold_poisson_gaussian,dVdt}; + + const SmartCreateFactory create_factory_ele(*pc_product_ele); + const SmartCreateFactory create_factory_pos(*pc_product_pos); + const auto CreateEle = create_factory_ele.getSmartCreate(); + const auto CreatePos = create_factory_pos.getSmartCreate(); + + const auto Transform = SchwingerTransformFunc{m_qed_schwinger_y_size, + ParticleStringNames::to_index.find("w")->second}; + + const auto num_added = filterCreateTransformFromFAB<1>( dst_ele_tile, + dst_pos_tile, box, array_EMFAB, np_ele_dst, + np_pos_dst,Filter, CreateEle, CreatePos, + Transform); + + setNewParticleIDs(dst_ele_tile, np_ele_dst, num_added); + setNewParticleIDs(dst_pos_tile, np_pos_dst, num_added); + + } +} + void MultiParticleContainer::doQedEvents() { WARPX_PROFILE("MPC::doQedEvents"); @@ -1151,6 +1286,17 @@ void MultiParticleContainer::CheckQEDProductSpecies() "ERROR: Quantum Synchrotron product species is of wrong type"); } } + + if (m_do_qed_schwinger) { + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + allcontainers[m_qed_schwinger_ele_product]-> + AmIA() + && + allcontainers[m_qed_schwinger_pos_product]-> + AmIA(), + "ERROR: Schwinger process product species are of wrong type"); + } + } #endif -- cgit v1.2.3