[Numpy-discussion] Re: [EXTERNAL] Copy on __setitem__

[Benjamin Root via NumPy-Discussion]

Seems to make sense to me? Or is the following a bug?

>>> import numpy as np
>>> u = np.zeros(5)
>>> v = np.ones(5)
>>> u
array([0., 0., 0., 0., 0.])
>>> u[...] = v
>>> u
array([1., 1., 1., 1., 1.])
>>> v[4] = 5
>>> v
array([1., 1., 1., 1., 5.])
>>> u
array([1., 1., 1., 1., 1.])

If you don't do a copy, then it is a view, right? And so, should modifying
v[4] change u[4]? Relatedly, if these were object arrays, of mutable
objects, should mutating u[4] point to the exact same instance as v[4] and
mutating it in one array means that the other array "sees" the same changes?

Sorry if I'm missing the point, just want to make sure we don't
inadvertently change important implicit behavior.

Ben Root

On Thu, Dec 26, 2024 at 7:07 AM Mateusz Sokol  wrote:

> It looks like general assignment function `array_assign_subscript` calls
> `PyArray_CopyObject(view, op)` link
> 
> which in turn always calls `PyArray_DiscoverDTypeAndShape` with `copy=1`
> and this is where `copy=True` takes place.
>
> We could add a `copy` parameter to `PyArray_CopyObject` or set `copy=-1`
> (to `None`) for the `PyArray_DiscoverDTypeAndShape` call.
>
> Mateusz
>
> On Thu, Dec 26, 2024 at 10:19 AM Sebastian Berg <
> sebast...@sipsolutions.net> wrote:
>
>> That seems like a bug but not sure why it would happen.  It needs to
>> call `__array__`, but indeed of course not with `copy=True`.
>>
>> Would you open an issue on github?
>>
>> - Sebastian
>>
>>
>> On Thu, 2024-12-26 at 03:46 +, Israel, Daniel M via NumPy-
>> Discussion wrote:
>> > Sure.  I didn’t originally, because I thought it would require an
>> > entire custom array container, but the following trivial example
>> > actually shows the behavior in question:
>> >
>> > import numpy
>> >
>> > class MyThing(object):
>> > def __array__(self, dtype=None, copy=None):
>> > print(f"MyThing.__array__(dtype={dtype}, copy={copy})")
>> > return numpy.ones((5, 5))
>> >
>> > u = numpy.zeros((5, 5))
>> > v = MyThing()
>> >
>> > u[...] = v
>> >
>> > If you run this code, as part of the final assignment statement, the
>> > __array__ method is called for ‘v’ with copy=True.  Why?
>> >
>> > —
>> > Daniel Israel
>> > XCP-4: Continuum Models and Numerical Algorithms
>> > d...@lanl.gov
>> >
>> > On Dec 25, 2024, at 3:23 PM, Steven Ellis
>> >  wrote:
>> >
>> > Hi David,
>> >
>> > New to the listserv, but, maybe you can provide a reproducible
>> > example?
>> >
>> > Steven
>> >
>> > On Wed, Dec 25, 2024, 2:19 PM Israel, Daniel M via NumPy-Discussion
>> > mailto:numpy-discussion@python.org>>
>> > wrote:
>> > I was updating some code that uses a custom array container built
>> > with the mixin library.  Specifically, I was trying to eliminate some
>> > warnings due to the change to the __array__ interface to add a copy
>> > argument.  In doing so, I discovered that, for two objects u, v in my
>> > container class, the code:
>> >
>> > u[…] = v
>> >
>> > performs a copy on v.  Specifically, it calls __array__() with
>> > copy=True.  This seems unnecessary and wasteful of memory.  Can
>> > someone explain to me what is happening here?
>> >
>> > Thanks.
>> >
>> > —
>> > Daniel Israel
>> > XCP-4: Continuum Models and Numerical Algorithms
>> > d...@lanl.gov
>> >
>> > ___
>> > NumPy-Discussion mailing list --
>> > numpy-discussion@python.org
>> > To unsubscribe send an email to
>> > numpy-discussion-le...@python.org> > numpy-discussion-le...@python.org>
>> > https://mail.python.org/mailman3/lists/numpy-discussion.python.org/<
>> >
>> https://urldefense.com/v3/__https://mail.python.org/mailman3/lists/numpy-discussion.python.org/__;!!Bt8fGhp8LhKGRg!A_N5Cns4800odNjjStN3QAH94TwUNBEkvr7lqW4-PnnUiGS_kckVyXKfS4QTR4YADAXXGOTbKZ-6RDYTC7knRCMj$
>> > >
>> > Member address:
>> > stevenalonzoel...@gmail.com
>> >
>> > ___
>> > NumPy-Discussion mailing list -- numpy-discussion@python.org
>> > To unsubscribe send an email to numpy-discussion-le...@python.org
>> > https://mail.python.org/mailman3/lists/numpy-discussion.python.org/
>> > Member address: sebast...@sipsolutions.net
>>
>> ___
>> NumPy-Discussion mailing list -- numpy-discussion@python.org
>> To unsubscribe send an email to numpy-discussion-le...@python.org
>> https://mail.python.org/mailman3/lists/numpy-discussion.python.org/
>> Member address: mso...@quansight.com
>>
> ___
> NumPy-Discussion mailing list -- numpy-discussion@python.org
> To unsubscribe send an email to numpy-discussion-le...@python.org
> https://mail.python.org/mailman3/lists/numpy-discussion.python.org/
> Member address: ben.v.r...@gmail.com
>
___
NumPy-D

[Numpy-discussion] Re: slice an array along specific axis

[Benjamin Root via NumPy-Discussion]

When I usually need to do something like that, I just construct a tuple of
slice() objects. No need to use swapaxes(). Or am I missing something?

On Sat, Feb 1, 2025 at 10:24 AM Michael Mullen 
wrote:

> Hello,
>
> I am writing a class handling NumPy arrays, and basing it off some
> computations I have done in C++ using the Eigen library. For tensors whose
> length are only known at runtime, the Eigen chip method is useful for
> slicing a tensor along a specific axis while keeping all other axes the
> same. I have not found a similar method in NumPy, but a simple
> implementation is
>
> def chip(A, axis, vals):
> return np.swapaxes(np.swapaxes(A, axis, -1)[..., vals], -1, axis)
>
> Example usage would be (to slice the 3rd axis of a 4D array):
> A = np.random.randn(10,11,12,13)
> B = chip(A, axis=2, vals=np.arange(0,6))
> B.shape #(10, 11, 6, 13)
>
> Since this may be useful for others, despite its simplicity, I thought it
> may be useful to have something similar in NumPy.
>
> Best,
> Mike
> ___
> NumPy-Discussion mailing list -- numpy-discussion@python.org
> To unsubscribe send an email to numpy-discussion-le...@python.org
> https://mail.python.org/mailman3/lists/numpy-discussion.python.org/
> Member address: ben.v.r...@gmail.com
>
___
NumPy-Discussion mailing list -- numpy-discussion@python.org
To unsubscribe send an email to numpy-discussion-le...@python.org
https://mail.python.org/mailman3/lists/numpy-discussion.python.org/
Member address: arch...@mail-archive.com

[Numpy-discussion] Re: Proposal to Extend NumPy Broadcasting for (D₁, D₂, ..., N, M) → (D₁, D₂, ..., K, M) When K is a Multiple of N (K % N == 0)

[Benjamin Root via NumPy-Discussion]

Shasang,

My main concern is that if there is a legitimate bug in someone's code that
is causing mismatched array sizes, this can mask that bug in certain
situations where the mismatch just so happens to produce arrays of certain
shapes. I am intrigued by the idea, though.

Ben Root

On Tue, Mar 25, 2025 at 9:29 AM Shasang Thummar via NumPy-Discussion <
numpy-discussion@python.org> wrote:

> Dear NumPy Developers,
>
>
> I hope you are doing well. I am writing to propose an enhancement to NumPy’s 
> broadcasting mechanism that could make the library even more powerful, 
> intuitive, and flexible while maintaining its memory efficiency.
>
> ## **Current Broadcasting Rule and Its Limitation**
> As per NumPy's current broadcasting rules, if two arrays have different 
> shapes, the smaller array can be expanded **only if one of its dimensions is 
> 1**. This allows memory-efficient expansion of data without copying. However, 
> if a dimension is greater than 1, NumPy **does not** allow expansion to a 
> larger size, even when the larger size is a multiple of the smaller size.
>
> ### **Example of Current Behavior (Allowed Expansion)**
> ```python
> import numpy as np
>
> A = np.array([[1, 2, 3]])  # Shape (1,3)
> B = np.array([[4, 5, 6],# Shape (2,3)
>   [7, 8, 9]])
>
> C = A + B  # Broadcasting works because (1,3) can expand to (2,3)
> print(C)
>
> *Output:*
>
> [[ 5  7  9]
>  [ 8 10 12]]
>
> Here, A has shape (1,3), which is automatically expanded to (2,3) without
> copying because *a dimension of size 1 can be stretched*.
>
> However, *NumPy fails when trying to expand a dimension where N > 1, even
> if the larger size is a multiple of N*.
> *Example of a Case That Fails (Even Though It Could Work)*
>
> A = np.array([[1, 2, 3],# Shape (2,3)
>   [4, 5, 6]])
>
> B = np.array([[10, 20, 30],  # Shape (4,3)
>   [40, 50, 60],
>   [70, 80, 90],
>   [100, 110, 120]])
>
> C = A + B  # Error! NumPy does not allow (2,3) to (4,3)
>
> *This fails with the error:*
>
> ValueError: operands could not be broadcast together with shapes (2,3) (4,3)
>
> *Why Should This Be Allowed?*
>
> If *a larger dimension is an exact multiple of a smaller one*, then
> logically, the smaller array can be *repeated* along that axis *without
> physically copying data*—just like NumPy does when broadcasting (1,M) to
> (N,M).
>
> In the above example,
>
>- A has shape (2,3), and B has shape (4,3).
>- Since 4 is *a multiple of* 2 (4 % 2 == 0), A could be *logically
>repeated 2 times* to become (4,3).
>- NumPy *already does* a similar expansion when a dimension is 1, so
>why not extend the same logic?
>
> *Proposed Behavior (Expanding N → M When M % N == 0)*
>
> Allow an axis with size N to be broadcast to size M *if and only if M is
> an exact multiple of N (M % N == 0)*. This is *just as memory-efficient*
> as the current broadcasting rules because it can be done using *stride
> tricks instead of copying data*.
> *Example of the Proposed Behavior*
>
> If NumPy allowed this new form of broadcasting:
>
> A = np.array([[1, 2, 3],# Shape (2,3)
>   [4, 5, 6]])
>
> B = np.array([[10, 20, 30],  # Shape (4,3)
>   [40, 50, 60],
>   [70, 80, 90],
>   [100, 110, 120]])
>
> # Proposed new broadcasting rule
> C = A + B
>
> print(C)
>
> *Expected Output:*
>
> [[ 11  22  33]
>  [ 44  55  66]
>  [ 71  82  93]
>  [104 115 126]]
>
> This works by *logically repeating A* to match B’s shape (4,3).
> *Why This is a Natural Extension of Broadcasting*
>
>- *Memory Efficiency:* Just like broadcasting (1,M) to (N,M), this
>expansion does *not* require physically copying data. Instead, strides
>can be adjusted to *logically repeat elements*, making this as
>efficient as current broadcasting.
>- *Intuitiveness:* Right now, broadcasting already surprises new
>users. If (1,3) can become (2,3), why not (2,3) to (4,3) when 4 is a
>multiple of 2? This would be a more intuitive rule.
>- *Extends Current Functionality:* This is *not* a new concept—it *extends
>the existing rule* where 1 can be stretched to any value. We are
>generalizing it to *any factor relationship* (N → M when M % N == 0).
>
> *Implementation Considerations*
>
> The logic behind NumPy’s current broadcasting already uses *stride tricks*
> for memory-efficient expansion. Extending it to handle (N, M) → (K, M)
> (where K % N == 0) would require:
>
>- Updating np.broadcast_shapes(), np.broadcast_to(), and related
>functions.
>- Extending the existing logic that handles expanding 1 to support
>factors as well.
>- Ensuring backward compatibility and maintaining performance.
>
> *Conclusion*
>
> I strongly believe this enhancement would make NumPy’s broadcasting *more
> powerful, more intuitive, and just as efficient as before*. The fact that
> we can implement this *without unnecessary copying* makes it a natural
> ex

[Numpy-discussion] Tricky ufunc implementation question

[Benjamin Root via NumPy-Discussion]

I'm looking at a situation where I like to wrap a C++ function that takes
two doubles as inputs, and returns an error code, a position vector, and a
velocity vector so that I essentially would have a function signature of
(N), (N) -> (N), (N, 3), (N, 3). When I try to use np.vectorize() or
np.frompyfunc() on the python version of this function, I keep running into
issues where it wants to make the outputs into object arrays of tuples. And
looking at utilizing PyUFunc_FromFuncAndData, it isn't clear to me how I
can tell it to expect those two output arrays to have a size 3 outer
dimension.

Are ufuncs the wrong thing here? How should I go about this? Is it even
possible?

Thanks in advance,
Ben Root
___
NumPy-Discussion mailing list -- numpy-discussion@python.org
To unsubscribe send an email to numpy-discussion-le...@python.org
https://mail.python.org/mailman3//lists/numpy-discussion.python.org
Member address: arch...@mail-archive.com

[Numpy-discussion] Re: Tricky ufunc implementation question

[Benjamin Root via NumPy-Discussion]

Warren,

I'm fine with implementing it in C. I just didn't think gufuncs were for
me. I couldn't tell from the description if it would be for my usecase
since I wasn't looping over subarrays, and I didn't see any good examples.
Maybe the documentation could be clearer. I'll have a look at your examples.

I did try that signature with np.vectorize() with the signature keyword
argument, but it didn't seem to work. Maybe it didn't work for the reasons
in that open issue.

Thank you,
Ben Root

On Fri, Jun 27, 2025 at 8:03 PM Warren Weckesser via NumPy-Discussion <
numpy-discussion@python.org> wrote:

> On Fri, Jun 27, 2025 at 5:29 PM Benjamin Root via NumPy-Discussion
>  wrote:
> >
> > I'm looking at a situation where I like to wrap a C++ function that
> takes two doubles as inputs, and returns an error code, a position vector,
> and a velocity vector so that I essentially would have a function signature
> of (N), (N) -> (N), (N, 3), (N, 3). When I try to use np.vectorize() or
> np.frompyfunc() on the python version of this function, I keep running into
> issues where it wants to make the outputs into object arrays of tuples. And
> looking at utilizing PyUFunc_FromFuncAndData, it isn't clear to me how I
> can tell it to expect those two output arrays to have a size 3 outer
> dimension.
> >
> > Are ufuncs the wrong thing here? How should I go about this? Is it even
> possible?
>
> Ben,
>
> It looks like the simplest signature for your core operation would be
> (),()->(),(3),(3), with broadcasting taking care of higher dimensional
> inputs.  Because not all the core shapes are scalars, that would
> require a *generalized* ufunc (gufunc).  There is an open issue
> (https://github.com/numpy/numpy/issues/14020) with a request for a
> function to generate a gufunc from a Python function.
>
> numba has the @guvectorize decorator, but I haven't use it much, and
> in my few quick attempts just now, it appeared to not accept fixed
> integer sizes in the output shape.  But wait to see if any numba gurus
> respond with a definitive answer about whether or not it can handle
> the shape signature (),()->(),(3),(3).
>
> You could implement the gufunc in a C or C++ extension module, if you
> don't mind the additional development effort and packaging hassle.  I
> know that works--I've implemented quite a few gufuncs in ufunclab
> (https://github.com/WarrenWeckesser/ufunclab).
>
> Warren
>
>
> >
> > Thanks in advance,
> > Ben Root
> > ___
> > NumPy-Discussion mailing list -- numpy-discussion@python.org
> > To unsubscribe send an email to numpy-discussion-le...@python.org
> > https://mail.python.org/mailman3//lists/numpy-discussion.python.org
> > Member address: warren.weckes...@gmail.com
> ___
> NumPy-Discussion mailing list -- numpy-discussion@python.org
> To unsubscribe send an email to numpy-discussion-le...@python.org
> https://mail.python.org/mailman3//lists/numpy-discussion.python.org
> Member address: ben.v.r...@gmail.com
>
___
NumPy-Discussion mailing list -- numpy-discussion@python.org
To unsubscribe send an email to numpy-discussion-le...@python.org
https://mail.python.org/mailman3//lists/numpy-discussion.python.org
Member address: arch...@mail-archive.com