See #15
Karen found a similar issue which went away after reverting to numpy 1.26 from numpy 2. numpy 2 does say it changes how dtypes are silently cast in the case of math involving different dtypes, so I need to track down the exact operation(s) where this is happening