Skip to content

Feature/scalar stats budgets#241

Merged
adperezm merged 15 commits intoExtremeFLOW:mainfrom
nick-morse:feature/scalar_stats_budgets
Feb 26, 2026
Merged

Feature/scalar stats budgets#241
adperezm merged 15 commits intoExtremeFLOW:mainfrom
nick-morse:feature/scalar_stats_budgets

Conversation

@nick-morse
Copy link
Contributor

Added scalar statistics and budget calculations from interpolated points.

The implementation is based on the RS_budgets routines for the corresponding fluid statistics/budgets. As such, the code is messy, and should be cleaned up at some point. But I think it's worth getting this in as a baseline before cleanup.

I verified that $\langle \theta \rangle$, $\langle \theta'^2 \rangle$, $\langle u_i' \theta' \rangle$, and budgets of $\langle \theta'^2 \rangle$, $\langle u' \theta' \rangle$, and $\langle v' \theta' \rangle$ from a channel flow at $Re_\tau = 180$, $Pr = 0.71$ averaged from a 2D spanwise plane match the results from Neko PR #2149, which were calculated using a 1D CSV scalar statistics file from Neko. Additionally, the mean scalar balance looks good as well.

Here are the results, with reference data in black lines, CSV processed data in colored lines, and current implementation as symbols. Note that a perfect match isn't expected to the reference data as discussed in the reference PR, this is just a test of the statistics processing.
image

image image image

Copy link

Copilot AI left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Pull request overview

Adds passive scalar statistics and derived budget terms, mirroring existing Reynolds-stress budget tooling and enabling postprocessing from interpolated point data.

Changes:

  • Introduces a non-MPI interpolated-points pipeline for passive scalar stats and scalar budget computation into an HDF5 “*_format” output.
  • Adds a passive-scalar counterpart to existing RS budget utilities: derivative helpers, 2D→3D stats conversion, additional-field generation, and point interpolation workflows.

Reviewed changes

Copilot reviewed 2 out of 2 changed files in this pull request and generated 4 comments.

File Description
pysemtools/postprocessing/statistics/passive_scalar_budgets_interpolatedPoints_notMPI.py Reads interpolated scalar stat fields, reshapes/averages/saves them, and computes scalar budgets using fluid stats.
pysemtools/postprocessing/statistics/passive_scalar_budgets.py Adds helper routines to compute/write additional scalar-stat fields, convert 2D stats to 3D, and interpolate scalar/stat fields onto points.

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.


# filling in the missing z velocity
z_vel = fld3d.registry['s37']
fld3d.add_field(comm, field_name="w", field=z_vel, dtype=np.single)
Copy link

Copilot AI Feb 24, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

convert_scalar_2Dstats_to_3D supports datatype='double', but the added w field is always registered as np.single, which silently downcasts double inputs and makes the output inconsistent with wdsz=8. Use data_type (computed from datatype) for the dtype argument here.

Suggested change
fld3d.add_field(comm, field_name="w", field=z_vel, dtype=np.single)
fld3d.add_field(comm, field_name="w", field=z_vel, dtype=data_type)

Copilot uses AI. Check for mistakes.

# Check if coordinates match
if not np.allclose(XYZ, XYZ_fluid):
sys.exit("EXYZ coordinates in scalar and fluid files do not match!")
Copy link

Copilot AI Feb 24, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Typo in the error message: "EXYZ coordinates..." looks unintended and makes the message harder to search/grepp. Should be "XYZ coordinates...".

Suggested change
sys.exit("EXYZ coordinates in scalar and fluid files do not match!")
sys.exit("XYZ coordinates in scalar and fluid files do not match!")

Copilot uses AI. Check for mistakes.
elif msh.gdim == 2:
scalar_deriv2.c1 = coef.dudxyz(scalar_deriv.c1, coef.drdx, coef.dsdx, coef.dtdx)
scalar_deriv2.c2 = coef.dudxyz(scalar_deriv.c2, coef.drdy, coef.dsdy, coef.dtdy)
scalar_deriv2.c3 = 0.0 * scalar_deriv2.c3
Copy link

Copilot AI Feb 24, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In the 2D branch, scalar_deriv2.c3 is referenced before it is ever assigned (0.0 * scalar_deriv2.c3), which will raise AttributeError when msh.gdim == 2. Initialize c3 from an existing array (e.g., 0.0 * scalar_deriv2.c2) or explicitly allocate it based on c1/c2 shapes, similar to compute_scalar_first_derivative.

Suggested change
scalar_deriv2.c3 = 0.0 * scalar_deriv2.c3
scalar_deriv2.c3 = 0.0 * scalar_deriv2.c2

Copilot uses AI. Check for mistakes.
Comment on lines +1146 to +1149
if this_ext != this_ext_check:
warnings.warn(
"File index of fname_stat and fname_mean differ! Hope you know what you are doing!"
)
Copy link

Copilot AI Feb 24, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warnings.warn(...) is used here, but the warnings module is never imported in this function scope, so this will raise NameError the first time the file index differs. Add import warnings (or avoid warnings.warn and use print consistently).

Copilot uses AI. Check for mistakes.
@Mr-2c
Copy link
Collaborator

Mr-2c commented Feb 25, 2026

one rather unimportant question: some of these functions seem to be duplicate from RS_budgets.py. I was just wondeering whether it makes sense to import them instead such that once we start changing things we need to change it only once. the whole RS_budgets should probably be rewritten at some point though, so maybe it is actually better to keep them separate such that changes in one does not break the other...

@Mr-2c
Copy link
Collaborator

Mr-2c commented Feb 25, 2026

the rest looks okay in my extra superficial reading and since the budgets seem to be correct it can be merged I think. maybe only one suggestion I have is to add a semi working example (for instance the code that was used to create the plots included with this PR along with the stat file) in the 'unindexed' folder of 'examples' for easy future reference, sharing, etc.

@adperezm
Copy link
Collaborator

the whole RS_budgets should probably be rewritten at some point though, so maybe it is actually better to keep them separate such that changes in one does not break the other...

Yes, we talked about this as well. I think for now we leave it as is. We will need to refactor those parts of the code soonish anyway.

the rest looks okay in my extra superficial reading and since the budgets seem to be correct it can be merged I think. maybe only one suggestion I have is to add a semi working example (for instance the code that was used to create the plots included with this PR along with the stat file) in the 'unindexed' folder of 'examples' for easy future reference, sharing, etc.

This is a good idea. If possible an example would be nice.

@nick-morse
Copy link
Contributor Author

the rest looks okay in my extra superficial reading and since the budgets seem to be correct it can be merged I think. maybe only one suggestion I have is to add a semi working example (for instance the code that was used to create the plots included with this PR along with the stat file) in the 'unindexed' folder of 'examples' for easy future reference, sharing, etc.

I agree it would be nice to have an example. What I used above is based on the turb_channel example in Neko, so it could be nice to update that with a scalar solve and statistics. But for the channel one would always output a 1D line instead of the 2D slice I processed here, so maybe it's best to upload the stats files (2D) and write a small postprocessing notebook.

@adperezm
Copy link
Collaborator

the rest looks okay in my extra superficial reading and since the budgets seem to be correct it can be merged I think. maybe only one suggestion I have is to add a semi working example (for instance the code that was used to create the plots included with this PR along with the stat file) in the 'unindexed' folder of 'examples' for easy future reference, sharing, etc.

I agree it would be nice to have an example. What I used above is based on the turb_channel example in Neko, so it could be nice to update that with a scalar solve and statistics. But for the channel one would always output a 1D line instead of the 2D slice I processed here, so maybe it's best to upload the stats files (2D) and write a small postprocessing notebook.

Yes, I am doing that now (the second option). After I add that, I will merge. For now I will not have the notebook, just the script. And then we fix it later.

@adperezm
Copy link
Collaborator

The example runs, and I added the data to the sem_data repo. So I am Happy merging this now. Thanks for the contribution and time checking the code to all involved :D

@adperezm adperezm merged commit 03ace15 into ExtremeFLOW:main Feb 26, 2026
10 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

5 participants