Repository navigation
Conversation
and pass through **kw_apply_ufunc
and other improvements (will put more detailed description in PR overview).
Sevans711
added
new feature
|
Check out this pull request on See visual diffs & provide feedback on Jupyter Notebooks. Powered by ReviewNB |
ASV BenchmarkingBenchmark Comparison ResultsBenchmarks that have improved:
Benchmarks that have stayed the same:
Benchmarks that have got worse:
|
|
Re: ASV benchmarks, I think it looks okay. I am not expecting any performance changes, and I'm seeing no changes for the most part. The only exceptions are a small (~20%) slowdown in the 480km gradient benchmark (but no significant slowdown in the 120km gradient benchmark, so I think this might just be noise?), and a small increase to track_peakmem_gradient (~10%) but that also seems to occur randomly on different PRs which aren't touching gradient at all, so I think it's okay to ignore here? Requesting reviews because things mostly look okay. The readthedocs build isn't building right now but it's not this PR's fault, readthedocs just had a random issue. If reviewers can review the non-readthedocs parts for now (i.e. everything except the vector_calculus.ipynb changes) and then check the readthedocs parts in a later review, that would be awesome! |
cmdupuis3
left a comment
cmdupuis3
left a comment
There was a problem hiding this comment.
I have two speed enhancements on the cmd/2d+calculus branch.
Other issues:
- Results now inherit the input's attrs via apply_ufunc/binary ops (keep_attrs on in current xarray), so gradient/curl/divergence outputs carry the input's long_name/standard_name, which mislabels the derivative.
- Unit issues:
divergence(units=...)is now ignored. Also, curl and divergence units are missing their parentheses e.g. m/s/m vs. (m/s)/m
also adds tests for what the "units" attr should contain.
See discussion in #1573 for details. Basically, DeprecationWarning is not visible to users by default, but FutureWarning is.
|
As a sidenote, netcdf writers seem to complain when attrs' values are booleans, so I'm not sure if {"divergence": True}, {"curl": True}, {"divergence": True} is a great convention for calculus operators' results' attrs... but that seems to be the convention on |
|
@Sevans711 Yeah, there could definitely be more issues there. CF-compliant attributes are a giant cluster, so even what you see here may be deceptively oversimplified.
|
netcdf doesn't like boolean-valued attrs.
To confirm, we're not trying to solve all of them as part of this PR, though, right? This feels like another PR's worth of work, or possibly multiple PRs if the goal is to make it apply across all operations which alter results' values throughout all of uxarray.
Fixed by updating naming scheme to
Fixed by adding a
Fixed, great catch, I previously suspected deep=False was the default, but that is wrong, the default is deep=True! |
@Sevans711 Yeah, we're not trying to solve all of them, but we should definitely avoid making it worse in the code this PR touches. |
cmdupuis3
left a comment
cmdupuis3
left a comment
There was a problem hiding this comment.
import uxarray as ux
psi = ux.tutorial.open_dataset("outCSne30-vortex")["psi"]
psi.uxgrid.sphere_radius = 6.371e6
# explicitly unnamed inputs
anon = psi.rename(None)
print(anon.curl(anon).name) # curl_None_None
print(anon.divergence(anon).name) # divergence_None_None
print(anon.scalardotgradient(anon, anon).name) # scalardotgradient_None_None_None
print(anon.gradient()["zonal_gradient"].attrs["long_name"]) # Zonal gradient of None
# the common way it happens: xarray drops the name when combining
# two differently named arrays
u, v = psi.rename("u"), psi.rename("v")
u_total = u + 2 * v.rename("v2")
print(u_total.name) # None
print(u_total.divergence(v).name) # divergence_None_v
print(u_total.divergence(v).attrs["long_name"]) # Divergence of (None, v)This snippet produces nonsense names. curl_None_None is produced on main, but this branch is causing that to happen with divergence and scalardotgradient.
| grad_units = f"{base_units}/m" if base_units else "1/m" | ||
| else: | ||
| grad_units = f"{base_units}/rad" if base_units else "1/rad" | ||
| base_units = data.attrs.get("units", "1") |
There was a problem hiding this comment.
| base_units = data.attrs.get("units", "1") | |
| base_units = str(data.attrs.get("units") or "").strip() or "1" |
There was a problem hiding this comment.
What does this fix? The .strip() is surprising to me and feels out of place here (I can't think of a reason a user would have "units" attr start/end with whitespace, so if they did, I would want to respect that, rather than strip it without understanding the reasoning). Also, I think the original is more readable. This feels like a Claude-ism to make the code unnecessarily more complicated…?
There was a problem hiding this comment.
Think about what the old code covers; if you have a unit attribute "", that was originally handled by the else, but now you get /m instead of 1/m.
There was a problem hiding this comment.
Ah, thank you for clarifying, I was missing the edge case of arr.attrs["units"]="" (or " " or " \t\n " or any other combinations of whitespace). I do agree in case of empty string or just whitespace, the result should be something like "1/m" not "/m". The main reason I missed this was the lack of regression tests and comments about it, so I added regression tests and a relevant comment.
However, I think the way the old code handled other edge cases like arr.attrs["units"]=False or arr.attrs["units"]=0 was extremely misleading; both of those would be converted to "1/m". I'm here for supporting the simple "units is an empty string" case as a natural extension of units being formed from a combination of strings, but I don't think uxarray should be silently converting units=False into units of 1/m. Since there were no regression tests or documentation about cases like that, my guess is that they were an accidental consequence rather than something intentional.
I also added regression tests corresponding to the cases of units attr being None, False, 0. In these cases I'm asserting the original object should be converted to a string, so the result becomes something like "None/m", "False/m", or "0/m" (none of which are particularly meaningful, but at least they are clearly wrong instead of silently assuming those values for units should be treated as equivalent to "1").
|
@cmdupuis3 Thank you for the review and for finding all of my small typos! Fixed all the typos as you suggested.
I'm not sure why you called them "nonsense"; None is a marker of something being missing and, as you say, this follows the pattern that currently exists on main for curl. Something like "curl_unnamedarray_unnamedarray" feels a bit too verbose to me for the In general I would guess users will often see name and long_name simultaneously when displaying the array. In that view you can easily discover that None in the name means the corresponding input had no name. For example: 
Aside from this, and one remaining point on the inline comment (#1823 (comment), re: automatically rechunking along the grid dimension), I believe everything has been addressed! |
|
@Sevans711 I think just saying |
True, I agree it's less cluttered to use "divergence", and that None isn't super informative here. How should we handle Since the names aren't consistent on |
|
@Sevans711 Yeah, it does look a bit complicated. |
Honestly, something like u.curl(v) really isn't my favorite syntax either… it's really not a python idiom as much as it is a uxarray idiom. I usually think of attributes similarly to ownership, so I still read this mentally as "u's curl of v", like "the curl of v belonging to u." Another perfectly valid Python syntax could have been something like But, of course, syntax changes like that aren't in scope of this PR and we probably wouldn't want to touch them at this point since the feature already exists! |
the "standard" way to name results in other packages seems to be either: - result.name = inputs' name if all inputs have same name, else None - result.name = name of operation which was performed, regardless of inputs' names. This change uses the second option since calculus operations significantly change the meaning of the values.
|
From internal discussion with @erogluorhan, it sounds like the best way forward here is for this PR to update curl's name to behave like divergence's name and scalardotgradient's name on My understanding is that the main deciding factor here was to consider what other software packages do. For operations on multiple arrays, other packages seem to either (1) use a name based solely on the operation performed without incorporating the inputs' names at all (matching the decision made above), or (2) use the same exact name as the inputs if all inputs have the same name, else None (but matching inputs' names doesn't make sense here because calculus operations significantly change the meaning of the values). I kept the inputs' name information in the result's attrs["long_name"] for now, because it feels nice to be able to see that info, and I don't think it hurts to include? E.g., if you're much later in a workflow and start to lose track of which variable is which, having informative attrs may help to avoid mistakes. With that, I think the only outstanding point is the lingering comment (#1823 (comment)) about the units attr, which I believe I've addressed above? Otherwise, I think this is ready for re-review (@cmdupuis3)! |
|
More units issues:
|
@cmdupuis3 thank you for your reply. For (1), I'm not sure what you mean. Do you mean that if the sphere_radius is 1.0, but scale_by_radius=True, then we're returning units of "{base_units}/m"? I think this matches the behavior on For (2), as I mentioned in the original comment thread above, #1823 (comment):
Can you give an example where (2) is clearly a bug? I'm pretty sure that currently, units like |
|
(1) Try this: import numpy as np
import uxarray as ux
uxgrid = ux.open_grid("test/meshfiles/mpas/QU/oQU480.231010.nc")
uxda = ux.UxDataArray(
np.zeros(uxgrid.n_face), dims=["n_face"], attrs={"units": "K"}, uxgrid=uxgrid
)
uxgrid.sphere_radius = 1.0
print(uxda.gradient()["zonal_gradient"].attrs["units"]) # K/m, expected K/rad(2) Preserving the whitespace in " s " is an issue. We should have one standard format for all outputs IMO. You can bet that without sanitizing these strings, there will eventually be 20 different formats mixed together in one file. These units can get propagated very far downstream, and I would imagine the unstripped whitespace could e.g. wreck a lot of downstream parsers if that's what they're doing. |
Closes #1723
Overview
Adds support for applying calculus operators to 2D+ arrays (more dims than just the grid dimension), by creating the new
_apply_1dfunc_with_grid_core_dim()helper method and using it inside_compute_gradient(). The operators which now support 2D+ arrays are gradient(), curl(), divergence(), and scalardotgradient().Tests note: the tests which prove that these changes are compatible with dask arrays remain commented out for now; they should be uncommented when #1798 merges. See #1823 (comment) for more details.
While working on this I made numerous miscellaneous improvements related to calculus methods, expanding the scope of this PR as follows.
Small expansions of scope (examples):
xr.set_options(display_expand_data=False)). The corresponding icon still shows/hides, the default state is just "hidden" instead of shown.~with\~in one markdown cell (starts with "Pure rotation is incompressible...") which was rendering locally with strikethrough for some text instead of tildes.Sidenote about the notebook: there are a few cells in the middle which I didn't touch but the diff shows they were edited; this was done by ruff formatting.
Small expansions of scope (functionality):
Tiny expansions of scope (functionality):
grad_u_zonalandgrad_v_meridionalinstead of the confusingu = grad_u_zonalandv = grad_v_meridional.Expected Usage
Exactly the same syntax as the operators have currently, but they can now be applied to arrays with more dimensions:
PR Checklist
General
Testing & Benchmarking
Documentation and Examples
docs/api.rst; internal (private) function names start with an underscore (_)AI Disclosure
AI Usage: Claude for discussions and some code suggestions; GitHub Copilot's inline code suggestions