Skip to content

Commit

Permalink
Merge pull request yt-project#4806 from neutrinoceros/numpy2/trapezoi…
Browse files Browse the repository at this point in the history
…d_renaming

RFC: replace vendored version of numpy.trapz with supported API in numpy 2
  • Loading branch information
matthewturk authored Feb 25, 2024
2 parents 4d2d372 + 9f36d8e commit f60342a
Showing 1 changed file with 8 additions and 97 deletions.
105 changes: 8 additions & 97 deletions yt/_maintenance/numpy2_compat.py
Original file line number Diff line number Diff line change
@@ -1,101 +1,12 @@
# vendor functions that were moved from numpy 1.x to scipy
# avoid deprecation warnings in numpy >= 2.0

import functools
from importlib.metadata import version

import numpy as np
from numpy.core import overrides
from packaging.version import Version

array_function_dispatch = functools.partial(
overrides.array_function_dispatch, module="yt"
)
NUMPY_VERSION = Version(version("numpy"))


def _trapezoid_dispatcher(y, x=None, dx=None, axis=None):
return (y, x)


# from numpy 1.25 (numpy.trapz), deprecated in numpy 2.0 the function is vendored to
# avoid adding a runtime dependency on scipy.integrate.trapezoid the name of the
# original numpy function is also avoided,
# see https://github.com/scipy/scipy/issues/12924
@array_function_dispatch(_trapezoid_dispatcher)
def trapezoid(y, x=None, dx=1.0, axis=-1):
r"""
Integrate along the given axis using the composite trapezoidal rule.
If `x` is provided, the integration happens in sequence along its
elements - they are not sorted.
Integrate `y` (`x`) along each 1d slice on the given axis, compute
:math:`\int y(x) dx`.
When `x` is specified, this integrates along the parametric curve,
computing :math:`\int_t y(t) dt =
\int_t y(t) \left.\frac{dx}{dt}\right|_{x=x(t)} dt`.
Parameters
----------
y : array_like
Input array to integrate.
x : array_like, optional
The sample points corresponding to the `y` values. If `x` is None,
the sample points are assumed to be evenly spaced `dx` apart. The
default is None.
dx : scalar, optional
The spacing between sample points when `x` is None. The default is 1.
axis : int, optional
The axis along which to integrate.
Returns
-------
trapezoid : float or ndarray
Definite integral of `y` = n-dimensional array as approximated along
a single axis by the trapezoidal rule. If `y` is a 1-dimensional array,
then the result is a float. If `n` is greater than 1, then the result
is an `n`-1 dimensional array.
See Also
--------
sum, cumsum
Notes
-----
Image [2]_ illustrates trapezoidal rule -- y-axis locations of points
will be taken from `y` array, by default x-axis distances between
points will be 1.0, alternatively they can be provided with `x` array
or with `dx` scalar. Return value will be equal to combined area under
the red lines.
References
----------
.. [1] Wikipedia page: https://en.wikipedia.org/wiki/Trapezoidal_rule
.. [2] Illustration image:
https://en.wikipedia.org/wiki/File:Composite_trapezoidal_rule_illustration.png
"""
y = np.asanyarray(y)
if x is None:
d = dx
else:
x = np.asanyarray(x)
if x.ndim == 1:
d = np.diff(x)
# reshape to correct shape
shape = [1] * y.ndim
shape[axis] = d.shape[0]
d = d.reshape(shape)
else:
d = np.diff(x, axis=axis)
nd = y.ndim
slice1 = [slice(None)] * nd
slice2 = [slice(None)] * nd
slice1[axis] = slice(1, None)
slice2[axis] = slice(None, -1)
try:
ret = (d * (y[tuple(slice1)] + y[tuple(slice2)]) / 2.0).sum(axis)
except ValueError:
# Operations didn't work, cast to ndarray
d = np.asarray(d)
y = np.asarray(y)
ret = np.add.reduce(d * (y[tuple(slice1)] + y[tuple(slice2)]) / 2.0, axis)
return ret
if NUMPY_VERSION >= Version("2.0.0dev0"):
from numpy import trapezoid as trapezoid # type: ignore [attr-defined]
else:
from numpy import trapz as trapezoid # noqa: F401

0 comments on commit f60342a

Please sign in to comment.