Numpy support

[andrewgsavage]

The goal is to have a scalar class, AffineScalarFunc/ufloat, and an array class, UArray. The UArray class should be used to store arrays of ufloats rather than an ndarray of ufloats as is currently the case. __array_ufunc__ and __array_function__ can be used to overide numpy functions, and can return ufloats or UArray, (or bool or whatever the result should be)

__array_ufunc__
https://numpy.org/doc/stable/reference/arrays.classes.html#numpy.class.__array_ufunc__
This lets uncertainties control how numpy ufuncs operate. I think __array_function__ takes precedence so it'll only be used for scalar ops once thats implemented.

UArray
This would be an array that holds the nominal values and linearcombinations. Might be easier to use an ndarray of ufloats to start with than to store the floats in an ndarray seperate to the linearcombinations.

__array_function__
https://numpy.org/doc/stable/reference/arrays.classes.html#numpy.class.__array_function__
This lets uncertainties control how most numpy functions operate, eg np.mean(UArray). I'm not sure whether it gets used in np.mean(arr_of_ufloats) The difficult part here is that for every numpy function we implement, we'll need to write a function that can deal with any shape of input arrays and kwargs. I think we've only got 1D function wrappers at present. I think a reasonable goal is to get the functions currently in unumpy working for 1D inputs, and be at a point where a user can submit a PR implementing or improving a function they need.

[newville]

The goal is to have a scalar class, AffineScalarFunc/ufloat, and an array class, UArray.

OK. It is not obvious to me that UFloat needs to work with Numpy, but it seems nice to have.

The UArray class should be used to store arrays of ufloats rather than an ndarray of ufloats as is currently the case.

Why is that? I'm not saying I object, I just don't see what that is "should". Does "arrays of ufloats rather than an ndarray of ufloats" mean "list of ufloats", or do you really mean arrays (https://docs.python.org/3/library/array.html)?

I could certainly be missing some use cases, but I think an important case to consider (and what I know I would probably "normally have" in my own work, and sort of how lmfit would deal with this) is that the user would have an ndarray of "nominal values" and also an ndarray of "standard errors", say from some measurement uncertainty or statistical analysis. So, having UArray be initialized with

 my_uarray = UArray(values_ndarray, stderr_ndarray)

(where the shape, size, and dtype of values_ndarray and stderr_ndarray are the same) seems like an important case to consider and support. And, just as UFloat has a float nominal_value and a float std_dev, UArray could have an ndarray (of dtype Float64/Float32/Complex128) for nominal_value and an ndarray (same shape and dtype) for std_dev.

That could make UFloat and UArray act very much like one another. And the derivatives/linear components are then using the same kind of "extended representation". "float" -> ndarray. I would guess that most of the internal math also then "just works".

The idea of using lists or arrays (but not ndarray) of ufloats seems pretty different from that. What am I missing?

[andrewgsavage]

By "arrays of ufloats rather than an ndarray of ufloats", I meant that the UArray class should be used rather than the current unumpy.uarray, which gives an ndarray of ufloat objects.

Yes, you can store the nominal and std like that. Then you also need a numpy object dtype array to store the LinearCombination or a some kind of LinearCombinationArray, so that the LinearCombination can be used in operations. Continuing your example, the LinearCombination is needed to get this behaviour:

v = z - y
v.std_dev == x.std_dev
True

Looping will be necessary to allow z[0] = ufloat(1,1), which I'd intuitively expect to work. Maybe there's a way to shortcut and avoid looping in some cases, although that's somewhat more difficult to implement so I would leave for a while.

[newville]

By "arrays of ufloats rather than an ndarray of ufloats", I meant that the UArray

Sorry for being dense, I guess I don't know what UArray should do - I thought that is what we were discussing. It is entirely possible that I am missing something that is obvious to you.

Yes, you can store the nominal and std like that. Then you also need a numpy object dtype array to store the LinearCombination or a some kind of LinearCombinationArray, so that the LinearCombination can be used in operations.

Couldn't the linear combinations could just have ndarrays as values instead of floats. Why not?

Continuing your example, the LinearCombination is needed to get this behaviour:

v = z - y
v.std_dev == x.std_dev
True

Or, maybe that could be

v = z - y
v.std_dev == x.std_dev
array([ True,  True, True, ....])

? That would seem more natural to me....

Looping will be necessary to allow z[0] = ufloat(1,1)

Well, maybe

z.nominal_value[0] = 1
z.std_dev[0]  = 1

would be more like setting elements of arrays?

But wait, wasn't there a goal to make UFloat un-settable after creation? Why would UArray be settable after creation?

[andrewgsavage]

Sorry for being dense, I guess I don't know what UArray should do - I thought that is what we were discussing. It is entirely possible that I am missing something that is obvious to you.

After rereading I realise how unclear that statement was!

Couldn't the linear combinations could just have ndarrays as values instead of floats. Why not?

Could you show how you'd implement this? This is the bit where I'm not seeing an obvious, efficient way to implement it. Hence why I'm thinking an ndarray of AffineScalarFunc.

array([ True,  True, True, ....])

? That would seem more natural to me....

Yes, the point I was trying to make is it that to perform the adding or subtracting operations, it needs to use the LinearCombinations

Looping will be necessary to allow z[0] = ufloat(1,1)

Well, maybe

z.nominal_value[0] = 1
z.std_dev[0]  = 1

would be more like setting elements of arrays?

Again you're missing the LinearCombinations, which would need something like

z._linear_part[0] = ...

Otherwise it'd consider the previous values uncertainty in future operations

But wait, wasn't there a goal to make UFloat un-settable after creation? Why would UArray be settable after creation?

A float is an immutable object that gets stored in an ndarray and you can assign new floats to elements in the ndarray. I'd expect the same behaviour from a UArray.

[newville]

Could you show how you'd implement this?

Start with UFloat. Replace all "coerce to float" to "allow ndarrays", and checking that "std_dev" and "nominal_value" have the same dtype and shape. Maybe allow "std_dev" to be a scalar, upgrading it to ndarray with a constant value. There is then no need to come up with an efficient implementation.

v = z - y
v.std_dev == x.std_dev

would then just give

array([ True,  True, True, ....])

A float is an immutable object that gets stored in an ndarray and you can assign new floats to elements in the ndarray. I'd expect the same behaviour from a UArray.

Yes, it is just a question of what the representation of UArray is. You seem to insist that is correct and obvious that the only option is that UArray is some sequence of UFloat objects. Again, you have been remarkably vague in using "array" multiple times. In this discussion, the term "array" could have several very different technical interpretations. Do you mean Python list, Python array, on Numpy Array of object type? Yes, answer that question. Please do not ignore when someone asks you a direct question.

What I am saying is that UArray could be an object which contains a nominal value and standard deviation that are each numpy arrays of a numerical/float-like(int,float, maybe complex) dtype. That could make UArray be very much like UFloat in operations, where all math done to find and use derivatives and combine values to give Affine{Scalar/Array}Func objects are done with numpy functions (for, say sqrt) instead of "Python math".

And, in that case, yes numpy ndarrays are mutable one could assign to elements of them. You would do that as

z.nominal_value[0] = 1
z.std_dev[0]  = 1

[andrewgsavage]

Could you show how you'd implement this?

Start with UFloat. Replace all "coerce to float" to "allow ndarrays", and checking that "std_dev" and "nominal_value" have the same dtype and shape. Maybe allow "std_dev" to be a scalar, upgrading it to ndarray with a constant value. There is then no need to come up with an efficient implementation.

v = z - y
v.std_dev == x.std_dev

would then just give

array([ True,  True, True, ....])

That implementation will work for most operations, but there are some, like cumsum or setting an element of UArray to a new UFloat (ie with a linear combo that is independent to the other values in the array) where it won't behave as I'd expect.

Yes, it is just a question of what the representation of UArray is. You seem to insist that is correct and obvious that the only option is that UArray is some sequence of UFloat objects. Again, you have been remarkably vague in using "array" multiple times. In this discussion, the term "array" could have several very different technical interpretations. Do you mean Python list, Python array, on Numpy Array of object type? Yes, answer that question. Please do not ignore when someone asks you a direct question.

I am using array or ndarray as a shorthand for a numpy-like array that works with numpy. It is vague because it is not yet defined how they will be stored. I have not mentioned lists or the array module.

What I am saying is that UArray could be an object which contains a nominal value and standard deviation that are each numpy arrays of a numerical/float-like(int,float, maybe complex) dtype. That could make UArray be very much like UFloat in operations, where all math done to find and use derivatives and combine values to give Affine{Scalar/Array}Func objects are done with numpy functions (for, say sqrt) instead of "Python math".

And, in that case, yes numpy ndarrays are mutable one could assign to elements of them. You would do that as

z.nominal_value[0] = 1
z.std_dev[0]  = 1

This is where I got hung up on the implementation so went to the numpy array of objects. z.std_dev[0] = 1 Here you have a UArray and are setting its first value. Does this mean you are expanding the UArray's linear combination to create a numpy array of floats for its std_dev and then setting the first value? Then how do you deal with the std_dev in subsequent operations? By setting the std_dev for the first element, it's linear combination should be different to other elements.

Nevertheless your implemention does should work for most operations. Perhaps we should have two UArray classes then, one more perfomant but with some unsupported capabilities.

[jagerber48]

We're all aware that this code works:

import numpy as np          
from uncertainties import uf
                            
x = ufloat(10, 1)           
y = ufloat(20, 2)           
z = ufloat(30, 3)           
                            
arr_1 = np.array([x, y, z]) 
print(f'{arr_1=}')          
                            
arr_2 = 2 * arr_1           
print(f'{arr_2-arr_1-arr_1=}
print(np.mean(arr_1))       

# arr_1=array([10.0+/-1.0, 20.0+/-2.0, 30.0+/-3.0], dtype=object)
# arr_2-arr_1-arr_1=array([0.0+/-0, 0.0+/-0, 0.0+/-0], dtype=object)
# 20.0+/-1.2

This didn't require any dependency of uncertainties on numpy.

This code currently breaks

x = ufloat(10, 1)
print(f'{np.sqrt(x)=}')
# TypeError: loop of ufunc does not support argument 0 of type Variable which has no callable sqrt method

But if we do

Variable.sqrt = lambda x: ufloat(12345, 42)
x = ufloat(10, 1)
print(np.sqrt(x))

# (1.234+/-0.004)e+04

It works. Again we get numpy interoperability without actually needing to import numpy into uncertainties. We simply need to define sqrt on Variable or AffineScalarFunc. Presumably we could do similar on all the ufuncs defined on this page https://numpy.org/doc/stable/reference/ufuncs.html#available-ufuncs. In that way we could get numpy interoperability (the ability to create and perform mathematical operations on arrays of UFloat objects) without needing to actually import numpy into uncertainties. Perhaps using these ufunc hooks is a way to shortcut this process and essentially wrap all of the ready-built numpy ufuncs with one with a UFloat compatibility wrapper?? @andrewgsavage is this roughly the idea?

Note that defining sqrt on Variable doesn't let math.sqrt to work on UFloat since math.sqrt tries to cast its inputs to floats.

---

This would be very nice to get the ability to work with numpy arrays of UFloat without uncertainties even needing to import numpy!!

The only thing that is missing is the ability to conveniently construct arrays of UFloat from arrays of values and std_dev (or values and uncertainty linear combinations?) and the ability to conveniently deconstruct an array of UFloat into values and std_devs. I somehow doubt we can get this without import numpy into uncertainties. Perhaps there's a way to add a switch case in __array_ufunc__ that could handle functions we want to be added in when UFloat is turned into a numpy array... not sure.

[newville]

@jagerber48

Yeah, UArray as "ndarray of object type containing UFloats" sort of works, but I think that having a UArray that was extension of UFloat with ndarrays as nominal value and standard deviation would be both more useable and easier to maintain. Like, all of the standard python operators are already overwritten by Numpy -- we probably just need to determine whether we're using math.sqrt or numpy.sqrt , and so forth.

It might even be that they could be the same class UVariable, that can have nominal/uncertainty values that are a) coerced to be floats but can be (consistent) ndarrays if numpy is available.

And, yeah, we might need to check the math if extending to complex....

[jagerber48]

I would really appreciate seeing some test cases to get us on the same page about what the supported API should look like. I wonder if some of the confusion in this discussion is coming from differing assumptions about what API we’re aiming for.

I’ll try to write some tests this week to share or would look at what others can write in the mean time.

[newville]

What about

• Constructing a UArray from an ndarray of UFloat?

That would be acceptable to me, but I think it would be a lot less common. How often does someone have a collection of scalar ufloats, and want to turn those into an ordered ndarray?

And, that ndarray has dtype of "object. I'll admit that I have a dim view numpy "object" dtype: "object" ndarrays report to be an ndarray, and you can manipulate the shapes and do fancy indexing/slicing, but many operations that work for numerical ndarrarys fail.

• Doing ndarray - UArray or scalar - UArray operations?

That should be fine. If the nominal values and uncertainties are "floats or ndarrays", I think that will just work with extra steps, just as "scalar + ndarray" works in Python. The code for UFloat and UArray could be nearly identical, if not actually the same.

• Using numpy function on a UArray: e.g. np.sin(), np.mean().

I think that, as shown very nicely by @andrewgsavage, that can be made to also "just work", and I think that would also extend readily from floats to ndarrays of floats. I think an important advantage of that approach is that it could work for other numerical functions, say from scipy.special.

Anyway, sorry that I do not have time at the moment to work on the coding part of this. But I think your approach is really promising.

[jagerber48]

What about

• Constructing a UArray from an ndarray of UFloat?

That would be acceptable to me, but I think it would be a lot less common. How often does someone have a collection of scalar ufloats, and want to turn those into an ordered ndarray?

I'm imagining a workflow where I have a lot of pieces of data that I need to do a few manipulations to as I go. So, if I'm working on a single piece of data at a time, I might generate a UFloat early in the process so I can do the intermediate steps with error propagation and only at the end bundle into a UArray. I admit it's possible that this workflow could possibly be refactored to use UArray for more of the process. Even if this latter case is the better/more performant way to do the workflow there is at least one user whose first inclination for some problems would be to end up with a bunch of UFloat objects that I want to bundle into a UArray. But it would be simple. Something like:

import numpy as np

class UArray:
    ...

    @classmethod
    def from_ndarray_of_ufloat(cls, nd_arr):
        val_arr = np.vectorize(lambda: x.n)(nd_arr)
        std_arr = np.vectorize(lambda: x.s)(nd_arr)
        return cls(val_arr, std_arr)

The question is just do we include this convenience method or not?

---

• Doing ndarray - UArray or scalar - UArray operations?

That should be fine. If the nominal values and uncertainties are "floats or ndarrays", I think that will just work with extra steps, just as "scalar + ndarray" works in Python. The code for UFloat and UArray could be nearly identical, if not actually the same.

agreed.

---

I think that, as shown very nicely by @andrewgsavage, that can be made to also "just work", and I think that would also extend readily from floats to ndarrays of floats.

Ok, this comment inspired me to look more deeply into #249. @andrewgsavage There's something basic I'm still not understanding. That PR defines UFloatNumpy. Is UFloatNumpy to be thought of as a scalar? That's what it looks like to me looking at the code and the tests. But, in the current state of that PR, is UFloatNumpy to be also thought of as an array? It doesn't look like it, but it does implement __array_ufunc__ and __array_function__. Are those methods meant to be implemented on scalars or arrays? From some quick testing on that branch, it looks to me like UFloatNumpy is like a scalar that plays better when inserted in numpy arrays than UFloat. np.sqrt will work on UFloatNumpy, and it will also work on an ndarray of UFloatNumpy. It looks like np.mean works on a ndarray array of UFloatNumpy.

However, I still don't think there's a way to access the nominal values or std_devs of an array of UFloatNumpy using an instance method like .n or .nominal_value etc. It seems like this is the only missing functionality that I see. So we would still need a helper function (rather than instance attribute) like def get_nominal_value(u_array): return np.vectorize(lambda x: x.n)(u_array).

---

Understanding better @andrewgsavage PR #249, I think I can start to understand the discussion about the implementation.

For @andrewgsavage I think the question is how can we get something like u_arr.nominal_value or u_arr.nominal_values to work where we get an array of nominal values from an array of or UFloat? Otherwise I think all the cases discussed here are satisfied.

@newville I'm not totally sure what your take on UFloatNumpy from #249 is, but maybe you're taking issue with (1) UFloatNumpy being implemented as a new class which is an optional parent of UFloat and (2) UFloatNumpy isn't an array, but is rather a scalar, so #249 still lacks a UArray class. It rather relies on ndarray of UFloat which you seem uneasy with for reasons like "but many operations that work for numerical ndarrarys fail." but I think this is exactly what is addressed by #249. That is, #249 makes it so that if we have an ndarray of UFloat (inherited from UFloatNumpy) then all (hopefully) of the numerical operations we care about now succeed.

---

edit: Actually #249 also lacks the ability to construct an array of UFloat/UFloatNumpy from two ndarray of float.

[newville]

@jagerber48 @andrewgsavage For sure, I am trying to take this opportunity as a fresh look to ask "what do we want uncertainties to do to support ndarrays"?

Sort of out-of-order:

For #249, I love that np.sin(ufloat) works. It's OK that math.sin(ufloat) cannot, and I get the np.sin() working sort of means that "hey maybe we don't need umath", at least if numpy is available.

But, if we have to use loops to support ndarrays, I think "there must be a better way".

I think we should not require numpy. But I also think that it is okay to say that a lot of stuff here does assume numpy. It is hard to imagine someone wanting to propagate uncertainties and not being able or willing to install numpy. So "if you want us to help you calculate the cosine of a value with uncertainties, install numpy" is not unreasonable.

For "object-ndarray of ufloats vs separate float-like (might be ints or complex) ndarrys for nominal and stderr", I think converting between such representations would be fine. Real float-like ndarrays are way better to use ;).

For the workflow example: OK, test/play around with Ufloats with float values, but then when you want to work with a real ndarray, it is almost certainly worth the effort to convert to float-like numpy ndarrays, and not some object thing that has to be iterated over.

[jagerber48]

Here's what I think we have for API requirements.

1. We want numpy operations like np.sin to work on a UFloat like scalar. I'll call this UScalar. It will probably be a parent, child, or equal to UFloat.
2. We want some kind of array of UScalars. I will call this UArray. It may be equal to UScalar, it might be a dtype=object ndarray of UScalar, or it might be some other kind of object that holds an array of float nominal values and uncertainty linear combinations. Note importantly that it does not suffice for the second array to be a float array of standard deviations. We need to keep track of correlations for each entry in the matrix.
3. It should be possible to construct a UArray from two dtype=float ndarrays: nominal value and std_dev arrays
4. It should be possible to construct a UArray from an array (ndarray, nested lists, etc.) of UScalar
5. It should be possible to do the following operations and gets results similar to ndarray
   a. UArray - UArray,
   b. UScalar-UArray,
   c. dtype=float ndarray - Uarray,
   e. float - Uarray
   f. I think numpy-like broadcasting etc. should be supported.
6. It should be possible to extract an array nominal_value, std_dev and the uncertainty linear combinations from a UArray.
7. numpy functions:
   a. Functions like np.sin which take an array and return a same-shaped array should work on UArray.
   b. Functions like np.mean which take an array and return a scalar should work on UArray.
   c. I also think numpy functions like np.mean(uarr, axis=1) which take an array and return a new array of a different shape should also work.

#249 largely covers point 1. Together with #249, Using dtype=object ndarray of UFloatNumpy (or UFloat) covers 2, 4, 5, and 7.

This leaves 3 and 6. What seems easiest to me is subclassing ndarray to make UArray and simply adding on methods to extract the nominal_value, std_dev and uncertainty linear combinations to cover 6 and setting up __init__ and some other user-accessible constructor methods to cover 2. I know there's some complexity that comes with subclassing ndarray, but it seems like we could overcome it and expose a nice API.

@newville do you agree with 1-7 as requirements?

---

Regarding changing UFloat to accept ndarrays as nominal_value and std_dev. I think there are challenges with this approach. First of all, passing std_dev alone doesn't give information about correlations, but that may be necessary. So we need to store the linear combination of uncertainty atoms (as they're called in my new branch, but the same concept exist in the current master branch) as well. This latter array with be dtype=object. I also don't see exactly what is gained by this approach of having UScalar = UArray. Thinking only about nominal values and std_dev and not thinking about correlations, and dropping a couple of the requirements above, I see how it could maybe be easily implemented. But getting it to hit all requirements may be tough.

[newville]

@jagerber48

1. We want numpy operations like np.sin to work on a UFloat like scalar. I'll call this UScalar. It will probably be a parent, child, or equal to UFloat.

OK.

2. We want some kind of array of UScalars.

Speaking only for myself, I want an object that is like UFloat but that uses 2 Numpy ndarrays of a numeric dtype (typically Float64) for both nominal value and standard error. A sequence (Python array, list, numpy object array) of UScalars that must be iterated over is a lot less useful.

I will call this UArray. It may be equal to UScalar, it might be a dtype=object ndarray of UScalar, or it might be some other kind of object that holds an array of float nominal values and uncertainty linear combinations.

OK. Again, to me, one of these is much more valuable.

Note importantly that it does not suffice for the second array to be a float array of standard deviations. We need to keep track of correlations for each entry in the matrix.

I think correlations between two UValues do not need to be stored differently in an individual UValue - that would work exactly as for UFloat, I think.

Non-trivial (ie, other than 1 or 0) correlations between two UValues could be used, if known, and stored for calculations with correlated UValues (float or ndarrays). I think this also suggests that using a "context" for calculations might be worth considering. Currently, context is only global (process-wide). Your approach using uuids for "Atoms" would possibly be "more global"..... But we could aim for a context in which a selected set of UValues exist and might have a correlation matrix....

But, correlations between different elements of a single UValue with an ndarray of values/uncertainties could be held in a separate container (say, "correl"). I don't see this as a terrible problem. I think it could also fall under "we will plan on supporting that eventually".

3.It should be possible to construct a UArray from two dtype=float ndarrays: nominal value and std_dev arrays
4. It should be possible to construct a UArray from an array (ndarray, nested lists, etc.) of UScalar

I agree. I do think it matters how they are stored and used internally. Using numerical ndarrays, avoiding loops will be very noticeable.

5. It should be possible to do the following operations and gets results similar to ndarray
   a. UArray - UArray,
   b. UScalar-UArray,
   c. dtype=float ndarray - Uarray,
   e. float - Uarray
   f. I think numpy-like broadcasting etc. should be supported.
6. It should be possible to extract an array nominal_value, std_dev and the uncertainty linear combinations from a UArray.

I agree with all of that.

7. numpy functions:
   a. Functions like np.sin which take an array and return a same-shaped array should work on UArray.
   b. Functions like np.mean which take an array and return a scalar should work on UArray.
   c. I also think numpy functions like np.mean(uarr, axis=1) which take an array and return a new array of a different shape should also work.

Yes. I think "a" (ufuncs that generate an ndarray of the same shape) should be straightforward (maybe not "easy").
I think "b" and "c" (reducing ufuncs) should not be too hard, though they might need some support code, especially with a view to handling element-to-element correlations.

Thanks!

[jagerber48]

Ok, I got what I consider to be a pretty easy implementation working on my branch: jagerber48#2

Here's the UArray class which is a subclass of ndarray. In this implementation a UArray is a dtype=object array of UFloats.
https://github.com/jagerber48/uncertainties/pull/2/files#diff-ed7b097e9a26f0dcd3627d33fffa6e9bd0581709bd4e1b987c5fb3380f93132bR11
It is a thin wrapper for ndarray that exposes nominal_value, std_dev, uncertainty properties and has from_val_arr_std_dev_arr alternative constructor to form a new array from val/std_dev arrays.

In umath_new.py I map many numpy ufunc names to umath functions and add these as attributes on the UFloat class.
https://github.com/jagerber48/uncertainties/pull/2/files#diff-0863db56ce7ba445f33177341e136f1b100a29f4a7e3cdc1e8002568382acd4fR124
This is what allows all the numpy ufuncs to naturally work on UFloat. Note that since __add__, sub, __mul__ etc. are already defined for UFloat it seems like the numpy ufuncs like np.add already defer to __add__.

Here's some tests that show we satisfy 1-7. I'd say in my opinion these tests pretty satisfactorily satisfy all of 1-7 except 7b which needs an asterisk.

import numpy as np

from uncertainties.core_new import UFloat
from uncertainties.unumpy.uarray import UArray

UFloat.__repr__ = lambda self: f'{self.n:.3f} +/- {self.s:.3f}'

x = UFloat(1, 0.1)
y = UFloat(2, 0.2)
z = UFloat(3, 0.3)

print("## Numpy operations (ufuncs) work on scalars (UFloat) ##")
print(f'{x=}')
print(f'{np.sin(x)=}')
print(f'{np.exp(x)=}')
print('')

print("## Constructing a UArray from an \"array\" of UFloat and looking at its properties ##")
uarr = UArray([x, y, z])
print(f'{uarr=}')
print(f'{uarr.nominal_value=}')
print(f'{uarr.std_dev=}')
print(f'{uarr.uncertainty=}')
print('')

print("## Constructing a UArray from 2 \"arrays\" of floats and looking at its properties ##")
uarr = UArray.from_val_arr_std_dev_arr([1, 2, 3], [0.1, 0.2, 0.3])
print(f'{uarr=}')
print(f'{uarr.nominal_value=}')
print(f'{uarr.std_dev=}')
print(f'{uarr.uncertainty=}')
print('')

print("## Binary operations with varying types")
narr = np.array([10, 20, 30])

print("# UArray : UArray")
print(f'{(uarr + uarr)=}')
print(f'{(uarr - uarr)=}')

print("# UArray : ndarray")
print(f'{(uarr + narr)=}')
print("# ndarray : UArray")
print(f'{(narr - uarr)=}')

print("# UFloat: UArray #")
print(f"{x * uarr}")
print("# UArray: UFloat #")
print(f"{uarr * x}")

print("# float : UArray #")
print(f"{42 * uarr}")
print("# UArray: float #")
print(f"{uarr * 42}")
print('')

print('## Numpy broadcasting works ##')
uarr1 = UArray.from_val_arr_std_dev_arr([[1, 2, 3], [4, 5, 6], [7, 8, 9]], np.ones((3, 3)))
uarr2 = UArray.from_val_arr_std_dev_arr([100, 1000, 1000], [10, 10, 10])
print(f'{uarr1=}')
print(f'{uarr2=}')
print(f'{(uarr1 + uarr2)=}')
print(f'{(uarr1 + uarr2).shape=}')
print('')

print('## More ufuncs work ##')

print('# np.mean #')
print(f'{np.mean(uarr1)=}')
print(f'{np.mean(uarr1, axis=0)=}')
print(f'{np.mean(uarr2)=}')

print('# other ufuncs #')
print(f'{np.exp(uarr1)=}')
print(f'{np.sin(uarr1)=}')
print(f'{np.sqrt(uarr1)=}')
print(f'{np.hypot(uarr1, uarr2)=}')
print(f'{np.hypot(uarr1, uarr2).shape=}')

Output:

## Numpy operations (ufuncs) work on scalars (UFloat) ##
x=1.000 +/- 0.100
np.sin(x)=0.841 +/- 0.054
np.exp(x)=2.718 +/- 0.272

## Constructing a UArray from an "array" of UFloat and looking at its properties ##
uarr=UArray([1.000 +/- 0.100, 2.000 +/- 0.200, 3.000 +/- 0.300], dtype=object)
uarr.nominal_value=array([1., 2., 3.])
uarr.std_dev=array([0.1, 0.2, 0.3])
uarr.uncertainty=array([((UncertaintyAtom(std_dev=0.1, uuid=UUID('cb411c52-c2ea-4c4a-99bd-3d21542108f4')), 1.0),),
       ((UncertaintyAtom(std_dev=0.2, uuid=UUID('4447b93f-c9e5-45d2-ae13-5e08aac70b64')), 1.0),),
       ((UncertaintyAtom(std_dev=0.3, uuid=UUID('a99b53f5-8fbc-4d15-8c46-2cd94fdd43b9')), 1.0),)],
      dtype=object)

## Constructing a UArray from 2 "arrays" of floats and looking at its properties ##
uarr=UArray([1.000 +/- 0.100, 2.000 +/- 0.200, 3.000 +/- 0.300], dtype=object)
uarr.nominal_value=array([1., 2., 3.])
uarr.std_dev=array([0.1, 0.2, 0.3])
uarr.uncertainty=array([((UncertaintyAtom(std_dev=0.1, uuid=UUID('57460d8d-3c1c-41f0-9882-56398ad5cadf')), 1.0),),
       ((UncertaintyAtom(std_dev=0.2, uuid=UUID('21e2afc3-4a86-4e04-b266-fc7724d4d038')), 1.0),),
       ((UncertaintyAtom(std_dev=0.3, uuid=UUID('0175f77e-e596-4bbf-95b6-df97778498da')), 1.0),)],
      dtype=object)

## Binary operations with varying types
# UArray : UArray
(uarr + uarr)=UArray([2.000 +/- 0.200, 4.000 +/- 0.400, 6.000 +/- 0.600], dtype=object)
(uarr - uarr)=UArray([0.000 +/- 0.000, 0.000 +/- 0.000, 0.000 +/- 0.000], dtype=object)
# UArray : ndarray
(uarr + narr)=UArray([11.000 +/- 0.100, 22.000 +/- 0.200, 33.000 +/- 0.300],
       dtype=object)
# ndarray : UArray
(narr - uarr)=UArray([9.000 +/- 0.100, 18.000 +/- 0.200, 27.000 +/- 0.300], dtype=object)
# UFloat: UArray #
UArray([1.000 +/- 0.141 2.000 +/- 0.283 3.000 +/- 0.424])
# UArray: UFloat #
UArray([1.000 +/- 0.141 2.000 +/- 0.283 3.000 +/- 0.424])
# float : UArray #
UArray([42.000 +/- 4.200 84.000 +/- 8.400 126.000 +/- 12.600])
# UArray: float #
UArray([42.000 +/- 4.200 84.000 +/- 8.400 126.000 +/- 12.600])

## Numpy broadcasting works ##
uarr1=UArray([[1.000 +/- 1.000, 2.000 +/- 1.000, 3.000 +/- 1.000],
        [4.000 +/- 1.000, 5.000 +/- 1.000, 6.000 +/- 1.000],
        [7.000 +/- 1.000, 8.000 +/- 1.000, 9.000 +/- 1.000]], dtype=object)
uarr2=UArray([100.000 +/- 10.000, 1000.000 +/- 10.000, 1000.000 +/- 10.000],
       dtype=object)
(uarr1 + uarr2)=UArray([[101.000 +/- 10.050, 1002.000 +/- 10.050, 1003.000 +/- 10.050],
        [104.000 +/- 10.050, 1005.000 +/- 10.050, 1006.000 +/- 10.050],
        [107.000 +/- 10.050, 1008.000 +/- 10.050, 1009.000 +/- 10.050]],
       dtype=object)
(uarr1 + uarr2).shape=(3, 3)

## More ufuncs work ##
# np.mean #
np.mean(uarr1)=5.000 +/- 0.333
np.mean(uarr1, axis=0)=UArray([4.000 +/- 0.577, 5.000 +/- 0.577, 6.000 +/- 0.577], dtype=object)
np.mean(uarr2)=700.000 +/- 5.774
# other ufuncs #
np.exp(uarr1)=UArray([[2.718 +/- 2.718, 7.389 +/- 7.389, 20.086 +/- 20.086],
        [54.598 +/- 54.598, 148.413 +/- 148.413, 403.429 +/- 403.429],
        [1096.633 +/- 1096.633, 2980.958 +/- 2980.958,
         8103.084 +/- 8103.084]], dtype=object)
np.sin(uarr1)=UArray([[0.841 +/- 0.540, 0.909 +/- 0.416, 0.141 +/- 0.990],
        [-0.757 +/- 0.654, -0.959 +/- 0.284, -0.279 +/- 0.960],
        [0.657 +/- 0.754, 0.989 +/- 0.146, 0.412 +/- 0.911]], dtype=object)
np.sqrt(uarr1)=UArray([[1.000 +/- 0.500, 1.414 +/- 0.354, 1.732 +/- 0.289],
        [2.000 +/- 0.250, 2.236 +/- 0.224, 2.449 +/- 0.204],
        [2.646 +/- 0.189, 2.828 +/- 0.177, 3.000 +/- 0.167]], dtype=object)
np.hypot(uarr1, uarr2)=UArray([[100.005 +/- 10.000, 1000.002 +/- 10.000, 1000.004 +/- 10.000],
        [100.080 +/- 9.992, 1000.012 +/- 10.000, 1000.018 +/- 10.000],
        [100.245 +/- 9.976, 1000.032 +/- 10.000, 1000.040 +/- 10.000]],
       dtype=object)
np.hypot(uarr1, uarr2).shape=(3, 3)

On 7b I had to make a "bespoke" hack to UArray.mean so that calling np.mean on UArray returned a scalar UFloat rather than a singleton 0-dimensional UArray of UFloat. https://github.com/jagerber48/uncertainties/pull/2/files#diff-ed7b097e9a26f0dcd3627d33fffa6e9bd0581709bd4e1b987c5fb3380f93132bR39

[jagerber48]

It's a little funny that we can do

import numpy as np

from uncertainties.core_new import UFloat
from uncertainties import umath_new

x = UFloat(np.pi/2, 0.3)

print(f'{x=}')
print(f'{x.sin()=}')
print(f'{x.cos()=}')

and get

# x=UFloat(1.5707963267948966, 0.3)
# x.sin()=UFloat(1.0, 1.8369701987210297e-17)
# x.cos()=UFloat(6.123233995736766e-17, 0.3)

but is that actually a problem? It's arguably a feature. It's maybe more weird for functions that take two arguments:

y = UFloat(np.pi/4, 0.1)
print(f'{x.hypot(y)=}')
# x.hypot(y)=UFloat(1.7562036827601817, 0.27202941017470883)

We could opt to only add in these functions onto UArray if numpy is imported but then UFloat has different behaviors depending on the environment. That seems Very Bad (tm). Another option would be to have a UArrayNumpy scalar class which inherits from UFloat but includes these methods. Then when we construct UArray we could cast UFloat to UFloatNumpy. But then when we pull things out of UArray we would get UFloatNumpy, not UFloat which would just be a bit annoying. Seems like added complexity for no gain.