diff --git a/pint/numpy_func.py b/pint/numpy_func.py
index 5994461..aa03e18 100644
--- a/pint/numpy_func.py
+++ b/pint/numpy_func.py
@@ -416,7 +416,6 @@ matching_input_copy_units_output_ufuncs = [
 copy_units_output_ufuncs = ["ldexp", "fmod", "mod", "remainder"]
 op_units_output_ufuncs = {
     "var": "square",
-    "prod": "size",
     "multiply": "mul",
     "true_divide": "div",
     "divide": "div",
diff --git a/pint/quantity.py b/pint/quantity.py
index 1f18166..f5453d2 100644
--- a/pint/quantity.py
+++ b/pint/quantity.py
@@ -22,6 +22,7 @@ from typing import List
 from packaging import version
 
 from .compat import (
+    HAS_NUMPY_ARRAY_FUNCTION,
     NUMPY_VER,
     _to_magnitude,
     babel_parse,
@@ -112,6 +113,18 @@ def check_implemented(f):
     return wrapped
 
 
+def method_wraps(numpy_func):
+    if isinstance(numpy_func, str):
+        numpy_func = getattr(np, numpy_func, None)
+
+    def wrapper(func):
+        func.__wrapped__ = numpy_func
+
+        return func
+
+    return wrapper
+
+
 @contextlib.contextmanager
 def printoptions(*args, **kwargs):
     """Numpy printoptions context manager released with version 1.15.0
@@ -1696,6 +1709,23 @@ class Quantity(PrettyIPython, SharedRegistryObject):
 
         return np.dot(self, b)
 
+    @method_wraps("prod")
+    def prod(self, *args, **kwargs):
+        """ Return the product of quantity elements over a given axis
+
+        Wraps np.prod().
+        """
+        # TODO: remove after support for 1.16 has been dropped
+        if not HAS_NUMPY_ARRAY_FUNCTION:
+            raise NotImplementedError(
+                "prod is only defined for"
+                " numpy == 1.16 with NUMPY_ARRAY_FUNCTION_PROTOCOL enabled"
+                f" or for numpy >= 1.17 ({np.__version__} is installed)."
+                " Please try setting the NUMPY_ARRAY_FUNCTION_PROTOCOL environment variable"
+                " or updating your numpy version."
+            )
+        return np.prod(self, *args, **kwargs)
+
     def __ito_if_needed(self, to_units):
         if self.unitless and to_units == "radian":
             return

synchronize Quantity.prod with np.prod(q)
Continuing from #867: #1087 fixed the calculated units for `prod` when the `axis` or `where` (not both) parameters are passed, but only for the ufunc and not for the method (also, right now this will calculate the units using input and output sizes, but instead we could just use `q.shape[axis]`).

I looked into implementing this but it seems this is somewhat complicated: the way it works right now, `__getattr__` will get the method from the magnitude and wrap it using `_numpy_method_wrap` which in turn will call `get_op_output_unit`. To make this work I think we'd need to move the units calculation of `_prod` into a new function that accepts the shape and the args (`axis` and `where`) and try to modify `get_op_output_unit` to accept functions in addition to names. @jthielen, thoughts?
diff --git a/pint/testsuite/test_numpy.py b/pint/testsuite/test_numpy.py
index 9ecdcb1..4a6b7e2 100644
--- a/pint/testsuite/test_numpy.py
+++ b/pint/testsuite/test_numpy.py
@@ -283,8 +283,14 @@ class TestNumpyMathematicalFunctions(TestNumpyMethods):
 
     # Sums, products, differences
 
+    @helpers.requires_array_function_protocol()
     def test_prod(self):
-        self.assertEqual(self.q.prod(), 24 * self.ureg.m ** 4)
+        axis = 0
+        where = [[True, False], [True, True]]
+
+        self.assertQuantityEqual(self.q.prod(), 24 * self.ureg.m ** 4)
+        self.assertQuantityEqual(self.q.prod(axis=axis), [3, 8] * self.ureg.m ** 2)
+        self.assertQuantityEqual(self.q.prod(where=where), 12 * self.ureg.m ** 3)
 
     @helpers.requires_array_function_protocol()
     def test_prod_numpy_func(self):
@@ -783,9 +789,6 @@ class TestNumpyUnclassified(TestNumpyMethods):
         self.assertQuantityAlmostEqual(np.std(self.q), 1.11803 * self.ureg.m, rtol=1e-5)
         self.assertRaises(OffsetUnitCalculusError, np.std, self.q_temperature)
 
-    def test_prod(self):
-        self.assertEqual(self.q.prod(), 24 * self.ureg.m ** 4)
-
     def test_cumprod(self):
         self.assertRaises(DimensionalityError, self.q.cumprod)
         self.assertQuantityEqual((self.q / self.ureg.m).cumprod(), [1, 2, 6, 24])

pip install -e .[dev]
parse_log_pytest
pytest
numpy>=1.14 pandas pytest pytest-mpl pytest-cov pytest-xdist pytest-mock pytest-asyncio
apt-get update apt-get install -y gcc
3.7
pytest --no-header -rA --tb=line --color=no -p no:cacheprovider -W ignore::DeprecationWarning
attrs @ file:///croot/attrs_1668696182826/work
certifi @ file:///croot/certifi_1671487769961/work/certifi
coverage==7.2.7
cycler==0.11.0
execnet==2.0.2
flit_core @ file:///opt/conda/conda-bld/flit-core_1644941570762/work/source/flit_core
fonttools==4.38.0
importlib-metadata @ file:///tmp/build/80754af9/importlib-metadata_1648562407465/work
importlib-resources==5.12.0
iniconfig @ file:///home/linux1/recipes/ci/iniconfig_1610983019677/work
Jinja2==3.1.6
kiwisolver==1.4.5
MarkupSafe==2.1.5
matplotlib==3.5.3
numpy==1.21.6
packaging @ file:///croot/packaging_1671697413597/work
pandas==1.3.5
Pillow==9.5.0
-e git+https://github.com/hgrecco/pint.git@23bb974664f624b8f7aa3b9f325a21cb6d049fdc#egg=Pint
pluggy @ file:///tmp/build/80754af9/pluggy_1648042572264/work
py @ file:///opt/conda/conda-bld/py_1644396412707/work
pyparsing==3.1.4
pytest==7.1.2
pytest-asyncio==0.21.2
pytest-cov==4.1.0
pytest-mock==3.11.1
pytest-mpl==0.17.0
pytest-xdist==3.5.0
python-dateutil==2.9.0.post0
pytz==2025.2
six==1.17.0
tomli @ file:///opt/conda/conda-bld/tomli_1657175507142/work
typing_extensions @ file:///croot/typing_extensions_1669924550328/work
zipp @ file:///croot/zipp_1672387121353/work

name: pint
channels:
  - defaults
  - https://repo.anaconda.com/pkgs/main
  - https://repo.anaconda.com/pkgs/r
  - conda-forge
dependencies:
  - _libgcc_mutex=0.1=main
  - _openmp_mutex=5.1=1_gnu
  - attrs=22.1.0=py37h06a4308_0
  - ca-certificates=2025.2.25=h06a4308_0
  - certifi=2022.12.7=py37h06a4308_0
  - flit-core=3.6.0=pyhd3eb1b0_0
  - importlib-metadata=4.11.3=py37h06a4308_0
  - importlib_metadata=4.11.3=hd3eb1b0_0
  - iniconfig=1.1.1=pyhd3eb1b0_0
  - ld_impl_linux-64=2.40=h12ee557_0
  - libffi=3.4.4=h6a678d5_1
  - libgcc-ng=11.2.0=h1234567_1
  - libgomp=11.2.0=h1234567_1
  - libstdcxx-ng=11.2.0=h1234567_1
  - ncurses=6.4=h6a678d5_0
  - openssl=1.1.1w=h7f8727e_0
  - packaging=22.0=py37h06a4308_0
  - pip=22.3.1=py37h06a4308_0
  - pluggy=1.0.0=py37h06a4308_1
  - py=1.11.0=pyhd3eb1b0_0
  - pytest=7.1.2=py37h06a4308_0
  - python=3.7.16=h7a1cb2a_0
  - readline=8.2=h5eee18b_0
  - setuptools=65.6.3=py37h06a4308_0
  - sqlite=3.45.3=h5eee18b_0
  - tk=8.6.14=h39e8969_0
  - tomli=2.0.1=py37h06a4308_0
  - typing_extensions=4.4.0=py37h06a4308_0
  - wheel=0.38.4=py37h06a4308_0
  - xz=5.6.4=h5eee18b_1
  - zipp=3.11.0=py37h06a4308_0
  - zlib=1.2.13=h5eee18b_1
  - pip:
      - coverage==7.2.7
      - cycler==0.11.0
      - execnet==2.0.2
      - fonttools==4.38.0
      - importlib-resources==5.12.0
      - jinja2==3.1.6
      - kiwisolver==1.4.5
      - markupsafe==2.1.5
      - matplotlib==3.5.3
      - numpy==1.21.6
      - pandas==1.3.5
      - pillow==9.5.0
      - pint==0.14.dev1+g23bb974
      - pyparsing==3.1.4
      - pytest-asyncio==0.21.2
      - pytest-cov==4.1.0
      - pytest-mock==3.11.1
      - pytest-mpl==0.17.0
      - pytest-xdist==3.5.0
      - python-dateutil==2.9.0.post0
      - pytz==2025.2
      - six==1.17.0
prefix: /opt/conda/envs/pint

pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_prod
pint/testsuite/test_numpy.py::TestNumpyArrayCreation::test_empty_like pint/testsuite/test_numpy.py::TestNumpyArrayCreation::test_full_like pint/testsuite/test_numpy.py::TestNumpyArrayCreation::test_ones_like pint/testsuite/test_numpy.py::TestNumpyArrayCreation::test_zeros_like pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_append pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_astype pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_atleast_1d pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_atleast_2d pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_atleast_3d pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_block_column_stack pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_broadcast_to pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_concat_stack pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_expand_dims pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_flat pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_flatten pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_flip_numpy_func pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_item pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_moveaxis pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_ravel pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_ravel_numpy_func pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_reshape pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_rollaxis pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_squeeze pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_swapaxes pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_transpose pint/testsuite/test_numpy.py::TestNumpyArrayManipulation::test_transpose_numpy_func pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_addition_with_incompatible_scalar pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_addition_with_scalar pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_cbrt pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_cross pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_cumprod pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_cumprod_numpy_func pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_diff pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_dot pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_dot_numpy_func pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_ediff1d pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_einsum pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_fix pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_gradient pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_nancumprod_numpy_func pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_nansum_numpy_func pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_power pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_prod_numpy_func pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_solve pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_sqrt pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_sum pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_sum_numpy_func pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_trapz pint/testsuite/test_numpy.py::TestNumpyMathematicalFunctions::test_unwrap pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_alen_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_all_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_allclose pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_any_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_argmax pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_argmax_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_argmin pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_argmin_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_argsort pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_argsort_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_array_protocol_unavailable pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_average_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_clip pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_clip_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_comparisons pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_compress pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_compress_nep18 pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_conj pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_copy_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_copyto pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_count_nonzero_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_cumprod pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_cumsum pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_cumsum_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_diagonal pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_diagonal_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_equal pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_fabs pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_fill pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_getitem pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_insert pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_interp_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_intersect1d pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_isclose_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_isin pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_iterable pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_iterator pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_max pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_max_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_max_with_axis_arg pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_max_with_initial_arg pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_maximum pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_mean pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_mean_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_median_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_meshgrid_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_min pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_min_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_min_with_axis_arg pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_min_with_initial_arg pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_minimum pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nan_to_num_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanargmax_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanargmin_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nancumsum_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanmax pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanmean_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanmedian_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanmin pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanpercentile pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanquantile pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanstd_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nanvar_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_ndarray_downcast pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_ndarray_downcast_with_dtype pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_ndim_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nonzero pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_nonzero_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_pad pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_percentile pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_pickle pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_ptp pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_ptp_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_put pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_quantile pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_repeat pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_resize pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_result_type_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_reversible_op pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_rot90 pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_round pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_round_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_searchsorted pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_searchsorted_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_setitem pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_shape pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_shape_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_sort pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_sort_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_std pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_std_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_take pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_tile pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_tolist pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_trace pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_trim_zeros_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_var pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_var_numpy_func pint/testsuite/test_numpy.py::TestNumpyUnclassified::test_where
pint/quantity.py pint/numpy_func.py
pint/quantity.py pint/numpy_func.py