-
Notifications
You must be signed in to change notification settings - Fork 201
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
[WIP] Add collisional diagnostic buffers #5543
base: development
Are you sure you want to change the base?
[WIP] Add collisional diagnostic buffers #5543
Conversation
This looks good. A comment - these new MultiFabs will automatically be accessible from Python since they are stored in the MultiFab map. You can also add convenience wrappers for them by adding the appropriate wrapper routines in |
for more information, see https://pre-commit.ci
@dpgrote, thanks for the tip! |
Source/Particles/Collision/BackgroundMCC/BackgroundMCCCollision.cpp
Outdated
Show resolved
Hide resolved
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Sorry for the slow response (again). Things have been quite hectic since the start of the year 😃
collision_energy_change, /**< Stores collisional energy transfer data for Monte Carlo Collisions */ | ||
collision_ionization /**< Stores ionization data for Monte Carlo Collisions */ |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
For clarity it would be good to add mcc
to these multifab names (to make it clear they are tied to the MCC routine). Maybe something like mcc_collision_energy_change
, etc.
if (WarpX::do_MCC_energy_tracking) { | ||
// Calculate the energy lost for tracking | ||
v_lost2 -= ux[ip]*ux[ip] + uy[ip]*uy[ip] + uz[ip]*uz[ip]; | ||
ParticleUtils::getEnergy(v_lost2, m, E_lost); | ||
|
||
// store the energy (add E_lost * w[ip] to the buffer of energy sent to the background) | ||
// get the particle's cell index (method used here is similar to TemperatureFunctor.cpp) | ||
const auto p = WarpXParticleContainer::ParticleType(ptd, ip); | ||
const auto [ii, jj, kk] = amrex::getParticleCell(p , plo, dxi).dim3(); | ||
|
||
// store energy lost in the buffer | ||
coll_E_change_arr(ii, jj, kk, 0) += E_lost * w[ip]; | ||
} |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It would be nice if this code does not need to be duplicated as it is currently. Perhaps the break
can be removed, and a else
used for the non-charge exchange processes?
// Control whether Monte Carlo Collision statistics are collected | ||
static bool do_MCC_energy_tracking; | ||
static bool do_MCC_ionization_tracking; | ||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
These should instead be made member variables of the BackgroundMCCCollision
class.
@@ -2224,6 +2274,25 @@ WarpX::AllocLevelMFs (int lev, const BoxArray& ba, const DistributionMapping& dm | |||
m_fields.alloc_init(FieldType::current_fp, Direction{1}, lev, amrex::convert(ba, jy_nodal_flag), dm, ncomps, ngJ, 0.0_rt); | |||
m_fields.alloc_init(FieldType::current_fp, Direction{2}, lev, amrex::convert(ba, jz_nodal_flag), dm, ncomps, ngJ, 0.0_rt); | |||
|
|||
// Allocate extra multifabs needed for MCC collision tracking | |||
if (do_MCC_ionization_tracking || do_MCC_energy_tracking) { |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
These multifabs could also be initialized in the constructor for BackgroundMCCCollision
.
Source/Particles/Collision/BackgroundMCC/BackgroundMCCCollision.cpp
Outdated
Show resolved
Hide resolved
Source/Particles/Collision/BackgroundMCC/BackgroundMCCCollision.cpp
Outdated
Show resolved
Hide resolved
def MCCEnergyWrapper(level=0, include_ghosts=False): | ||
return _MultiFABWrapper( | ||
mf_name="collision_energy_change", level=level, include_ghosts=include_ghosts | ||
) | ||
|
||
|
||
def MCCIonizWrapper(level=0, include_ghosts=False): | ||
return _MultiFABWrapper( | ||
mf_name="collision_ionization", level=level, include_ghosts=include_ghosts | ||
) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This may be a question for @dpgrote: Is this line in fields.py
all I need to do for multifab access in python? When I try the following in a picmi script:
ioniz_tracker = fields.MCCIonizWrapper()
tracking_data = ioniz_tracker[...]
I getting the following error:
File "/scratch/gpfs/user/warpx/Feb_2025/mcc_diag_check/Turner_inputs.py", line 425, in run_sim
tracking_data = ioniz_tracker[...]
File "/home/user/src/warpx/build_MCC_coll/lib/site-packages/pywarpx/fields.py", line 464, in __getitem__
result_global[global_slices] = f_arr
~~~~~~~~~~~~~^^^^^^^^^^^^^^^
ValueError: could not broadcast input array from shape (0,1,1,1) into shape (16,1,1,1)
Do you know how to fix this error? It seems like the multifab is either not being accessed correctly, or initialized correctly.
In some rudimentary checking, the multifabs seem to be correctly allocated (in this 1D case, 128 cells spread across 8 boxes of 16 for 8 MPI ranks) on the C++ side. It also looks like data is being saved in the multifab arrays correctly, as I can see the multifab array elements increase as ionization events occur (via print statements inserted in BackgroundMCCCollisions.cpp).
Stemming out of conversation #5381 with @roelof-groenewald, included is initial work on adding
multifab
buffers to track collisional diagnostics in a MCC-type simulation.The work on two buffers is in progress--one for ionization and one for total collisional energy loss. Both
multifabs
have been included inFields.H
and initialized inWarpX.cpp
. Relevant flags to turn on MCC tracking have been created inWarpX.H
and initialized inWarpX.cpp
. The energy lossmultifab
has been added into the collision protocol, and future commits will add the ionization buffer.In the end, I aim to set up Python access to these
multifabs
, with similar functionality to theparticle_containers.ParticleBoundaryBufferWrapper()
wrappers.