Feat: Adding Linear Algebra Dot operation support

quaddtype/README.md

@@ -25,11 +25,11 @@ np.array([1,2,3], dtype=QuadPrecDType("longdouble"))
 
 ## Installation from source
 
-The code needs the quad precision pieces of the sleef library, which
-is not available on most systems by default, so we have to generate
-that first. The below assumes one has the required pieces to build
-sleef (cmake and libmpfr-dev), and that one is in the package
-directory locally.
+The code needs the quad precision pieces of the sleef library, which is not available on most systems by default, so we have to generate that first. Choose the appropriate section below based on your operating system.
+
+### Linux/Unix/macOS
+
+The below assumes one has the required pieces to build sleef (cmake and libmpfr-dev), and that one is in the package directory locally.
 
 ```bash
 git clone --branch 3.8 https://github.com/shibatch/sleef.git
@@ -68,3 +68,91 @@ python -m pip install . -v --no-build-isolation -Cbuilddir=build -C'compile-args
 cd ..
 python -m pytest
 ```
+
+### Windows
+
+#### Prerequisites
+
+- **Visual Studio 2017 or later** (with MSVC compiler)
+- **CMake** (≥3.15)
+- **Python 3.10+** 
+- **Git**
+
+#### Step-by-Step Installation
+
+1. **Setup Development Environment**
+
+   Open a **Developer Command Prompt for VS** or **Developer PowerShell for VS** to ensure MSVC is properly configured.
+
+2. **Clone and Build SLEEF**
+
+   ```powershell
+   # Clone SLEEF library
+   git clone --branch 3.8 https://github.com/shibatch/sleef.git
+   cd sleef
+
+   # Configure with CMake for Windows
+   cmake -S . -B build -G "Visual Studio 17 2022" -A x64 -DSLEEF_BUILD_QUAD:BOOL=ON -DSLEEF_BUILD_SHARED_LIBS:BOOL=ON -DCMAKE_POSITION_INDEPENDENT_CODE=ON
+
+   # Build and install SLEEF
+   cmake --build build --config Release
+   cmake --install build --prefix "C:/sleef" --config Release
+
+   cd ..
+   ```
+
+3. **Setup Python Environment**
+
+   ```powershell
+   # Create and activate virtual environment
+   python -m venv numpy_quad_env
+   .\numpy_quad_env\Scripts\Activate.ps1
+
+   # Install build dependencies
+   pip install -U pip
+   pip install meson-python numpy pytest ninja meson
+   ```
+
+4. **Set Environment Variables**
+
+   ```powershell
+   # Set up paths and compiler flags
+   $env:INCLUDE = "C:/sleef/include;$env:INCLUDE"
+   $env:LIB = "C:/sleef/lib;$env:LIB"
+   $env:PATH = "C:/sleef/bin;$env:PATH"
+
+   # Note: QBLAS is disabled on Windows due to MSVC compatibility issues
+   $env:CFLAGS = "/IC:/sleef/include /DDISABLE_QUADBLAS"
+   $env:CXXFLAGS = "/IC:/sleef/include /DDISABLE_QUADBLAS"
+   $env:LDFLAGS = "C:/sleef/lib/sleef.lib C:/sleef/lib/sleefquad.lib"
+   ```
+
+5. **Build and Install numpy-quaddtype**
+
+   ```powershell
+   # Ensure submodules are initialized
+   git submodule update --init --recursive
+
+   # Build and install the package
+   python -m pip install . -v --no-build-isolation -Cbuilddir=build -C'compile-args=-v'
+   ```
+
+6. **Test Installation**
+
+   ```powershell
+   # Run tests
+   pytest -s tests/
+   ```
+
+1. **QBLAS Disabled**: QuadBLAS optimization is automatically disabled on Windows builds due to MSVC compatibility issues. This is handled by the `-DDISABLE_QUADBLAS` compiler flag.
+
+2. **Visual Studio Version**: The instructions assume Visual Studio 2022. For other versions, adjust the generator string:
+   - VS 2019: `"Visual Studio 16 2019"`
+   - VS 2017: `"Visual Studio 15 2017"`
+
+3. **Architecture**: The instructions are for x64. For x86 builds, change `-A x64` to `-A Win32`.
+
+4. **Alternative SLEEF Location**: If you prefer to install SLEEF elsewhere, update all path references accordingly.
+
+#### Windows Troubleshooting
+- **Link errors**: Verify that `sleef.lib` and `sleefquad.lib` exist in `C:/sleef/lib/`

quaddtype/numpy_quaddtype/__init__.py

@@ -8,8 +8,6 @@
     get_quadblas_version
 )
 
-import multiprocessing
-
 __all__ = [
     'QuadPrecision', 'QuadPrecDType', 'SleefQuadPrecision', 'LongDoubleQuadPrecision',
     'SleefQuadPrecDType', 'LongDoubleQuadPrecDType', 'is_longdouble_128', 

quaddtype/numpy_quaddtype/src/umath/matmul.cpp

@@ -36,7 +36,7 @@ quad_matmul_resolve_descriptors(PyObject *self, PyArray_DTypeMeta *const dtypes[
     if (descr_in1->backend != BACKEND_SLEEF || descr_in2->backend != BACKEND_SLEEF) {
         PyErr_SetString(PyExc_NotImplementedError,
                         "QBLAS-accelerated matmul only supports SLEEF backend. "
-                        "Other backends are not supported with QBLAS.");
+                        "Please raise the issue at SwayamInSync/QBLAS for longdouble support");
         return (NPY_CASTING)-1;
     }
 
@@ -61,7 +61,8 @@ quad_matmul_resolve_descriptors(PyObject *self, PyArray_DTypeMeta *const dtypes[
         QuadPrecDTypeObject *descr_out = (QuadPrecDTypeObject *)given_descrs[2];
         if (descr_out->backend != target_backend) {
             PyErr_SetString(PyExc_NotImplementedError,
-                            "QBLAS-accelerated matmul only supports SLEEF backend for output.");
+                        "QBLAS-accelerated matmul only supports SLEEF backend. "
+                        "Please raise the issue at SwayamInSync/QBLAS for longdouble support");
             return (NPY_CASTING)-1;
         }
         else {
@@ -118,7 +119,9 @@ quad_matmul_strided_loop_aligned(PyArrayMethod_Context *context, char *const dat
 
     QuadPrecDTypeObject *descr = (QuadPrecDTypeObject *)context->descriptors[0];
     if (descr->backend != BACKEND_SLEEF) {
-        PyErr_SetString(PyExc_RuntimeError, "Internal error: non-SLEEF backend in QBLAS matmul");
+        PyErr_SetString(PyExc_NotImplementedError,
+                        "QBLAS-accelerated matmul only supports SLEEF backend. "
+                        "Please raise the issue at SwayamInSync/QBLAS for longdouble support");
         return -1;
     }
 
@@ -206,7 +209,9 @@ quad_matmul_strided_loop_unaligned(PyArrayMethod_Context *context, char *const d
 
     QuadPrecDTypeObject *descr = (QuadPrecDTypeObject *)context->descriptors[0];
     if (descr->backend != BACKEND_SLEEF) {
-        PyErr_SetString(PyExc_RuntimeError, "Internal error: non-SLEEF backend in QBLAS matmul");
+        PyErr_SetString(PyExc_NotImplementedError,
+                        "QBLAS-accelerated matmul only supports SLEEF backend. "
+                        "Please raise the issue at SwayamInSync/QBLAS for longdouble support");
         return -1;
     }
 

quaddtype/tests/test_dot.py

@@ -1,109 +1,9 @@
 import pytest
 import numpy as np
+from test_utils import create_quad_array, assert_quad_equal, assert_quad_array_equal, arrays_equal_with_nan
 from numpy_quaddtype import QuadPrecision, QuadPrecDType
 
 
-# ================================================================================
-# UTILITIES
-# ================================================================================
-
-def assert_quad_equal(a, b, rtol=1e-15, atol=1e-15):
-    """Assert two quad precision values are equal within tolerance"""
-    # Ensure both operands are QuadPrecision objects for the comparison
-    if not isinstance(a, QuadPrecision):
-        a = QuadPrecision(str(a), backend='sleef')
-    if not isinstance(b, QuadPrecision):
-        b = QuadPrecision(str(b), backend='sleef')
-
-    # Use quad-precision arithmetic to calculate the difference
-    diff = abs(a - b)
-    tolerance = QuadPrecision(str(atol), backend='sleef') + QuadPrecision(str(rtol), backend='sleef') * max(abs(a), abs(b))
-
-    # Assert using quad-precision objects
-    assert diff <= tolerance, f"Values not equal: {a} != {b} (diff: {diff}, tol: {tolerance})"
-
-
-def assert_quad_array_equal(a, b, rtol=1e-25, atol=1e-25):
-    """Assert two quad precision arrays are equal within tolerance"""
-    assert a.shape == b.shape, f"Shapes don't match: {a.shape} vs {b.shape}"
-
-    flat_a = a.flatten()
-    flat_b = b.flatten()
-
-    for i, (val_a, val_b) in enumerate(zip(flat_a, flat_b)):
-        try:
-            assert_quad_equal(val_a, val_b, rtol, atol)
-        except AssertionError as e:
-            raise AssertionError(f"Arrays differ at index {i}: {e}")
-
-
-def create_quad_array(values, shape=None):
-    """Create a QuadPrecision array from values using Sleef backend"""
-    dtype = QuadPrecDType(backend='sleef')
-
-    if isinstance(values, (list, tuple)):
-        if shape is None:
-            # 1D array
-            quad_values = [QuadPrecision(str(float(v)), backend='sleef') for v in values]
-            return np.array(quad_values, dtype=dtype)
-        else:
-            # Reshape to specified shape
-            if len(shape) == 1:
-                quad_values = [QuadPrecision(str(float(v)), backend='sleef') for v in values]
-                return np.array(quad_values, dtype=dtype)
-            elif len(shape) == 2:
-                m, n = shape
-                assert len(values) == m * n, f"Values length {len(values)} doesn't match shape {shape}"
-                quad_matrix = []
-                for i in range(m):
-                    row = [QuadPrecision(str(float(values[i * n + j])), backend='sleef') for j in range(n)]
-                    quad_matrix.append(row)
-                return np.array(quad_matrix, dtype=dtype)
-
-    raise ValueError("Unsupported values or shape")
-
-
-def is_special_value(val):
-    """Check if a value is NaN or infinite"""
-    try:
-        float_val = float(val)
-        return np.isnan(float_val) or np.isinf(float_val)
-    except:
-        return False
-
-
-def arrays_equal_with_nan(a, b, rtol=1e-15, atol=1e-15):
-    """Compare arrays that may contain NaN values"""
-    if a.shape != b.shape:
-        return False
-
-    flat_a = a.flatten()
-    flat_b = b.flatten()
-
-    for i, (val_a, val_b) in enumerate(zip(flat_a, flat_b)):
-        # Handle NaN cases
-        if is_special_value(val_a) and is_special_value(val_b):
-            float_a = float(val_a)
-            float_b = float(val_b)
-            # Both NaN
-            if np.isnan(float_a) and np.isnan(float_b):
-                continue
-            # Both infinite with same sign
-            elif np.isinf(float_a) and np.isinf(float_b) and np.sign(float_a) == np.sign(float_b):
-                continue
-            else:
-                return False
-        elif is_special_value(val_a) or is_special_value(val_b):
-            return False
-        else:
-            try:
-                assert_quad_equal(val_a, val_b, rtol, atol)
-            except AssertionError:
-                return False
-
-    return True
-
-
 # ================================================================================
 # VECTOR-VECTOR DOT PRODUCT TESTS
 # ================================================================================
@@ -789,8 +689,4 @@ def test_dimension_mismatch_matrices(self):
         B = create_quad_array([1, 2, 3, 4, 5, 6], shape=(3, 2))  # Wrong size
 
         with pytest.raises(ValueError, match=r"matmul: Input operand 1 has a mismatch in its core dimension 0"):
-            np.matmul(A, B)
-
-
-if __name__ == "__main__":
-    pytest.main([__file__, "-v"])
+            np.matmul(A, B)

quaddtype/tests/test_utils.py

@@ -0,0 +1,98 @@
+from numpy_quaddtype import QuadPrecision, QuadPrecDType
+import numpy as np
+
+def assert_quad_equal(a, b, rtol=1e-15, atol=1e-15):
+    """Assert two quad precision values are equal within tolerance"""
+    # Ensure both operands are QuadPrecision objects for the comparison
+    if not isinstance(a, QuadPrecision):
+        a = QuadPrecision(str(a), backend='sleef')
+    if not isinstance(b, QuadPrecision):
+        b = QuadPrecision(str(b), backend='sleef')
+
+    # Use quad-precision arithmetic to calculate the difference
+    diff = abs(a - b)
+    tolerance = QuadPrecision(str(atol), backend='sleef') + QuadPrecision(str(rtol), backend='sleef') * max(abs(a), abs(b))
+
+    # Assert using quad-precision objects
+    assert diff <= tolerance, f"Values not equal: {a} != {b} (diff: {diff}, tol: {tolerance})"
+
+
+def assert_quad_array_equal(a, b, rtol=1e-25, atol=1e-25):
+    """Assert two quad precision arrays are equal within tolerance"""
+    assert a.shape == b.shape, f"Shapes don't match: {a.shape} vs {b.shape}"
+
+    flat_a = a.flatten()
+    flat_b = b.flatten()
+
+    for i, (val_a, val_b) in enumerate(zip(flat_a, flat_b)):
+        try:
+            assert_quad_equal(val_a, val_b, rtol, atol)
+        except AssertionError as e:
+            raise AssertionError(f"Arrays differ at index {i}: {e}")
+
+
+def create_quad_array(values, shape=None):
+    """Create a QuadPrecision array from values using Sleef backend"""
+    dtype = QuadPrecDType(backend='sleef')
+
+    if isinstance(values, (list, tuple)):
+        if shape is None:
+            # 1D array
+            quad_values = [QuadPrecision(str(float(v)), backend='sleef') for v in values]
+            return np.array(quad_values, dtype=dtype)
+        else:
+            # Reshape to specified shape
+            if len(shape) == 1:
+                quad_values = [QuadPrecision(str(float(v)), backend='sleef') for v in values]
+                return np.array(quad_values, dtype=dtype)
+            elif len(shape) == 2:
+                m, n = shape
+                assert len(values) == m * n, f"Values length {len(values)} doesn't match shape {shape}"
+                quad_matrix = []
+                for i in range(m):
+                    row = [QuadPrecision(str(float(values[i * n + j])), backend='sleef') for j in range(n)]
+                    quad_matrix.append(row)
+                return np.array(quad_matrix, dtype=dtype)
+
+    raise ValueError("Unsupported values or shape")
+
+
+def is_special_value(val):
+    """Check if a value is NaN or infinite"""
+    try:
+        float_val = float(val)
+        return np.isnan(float_val) or np.isinf(float_val)
+    except:
+        return False
+
+
+def arrays_equal_with_nan(a, b, rtol=1e-15, atol=1e-15):
+    """Compare arrays that may contain NaN values"""
+    if a.shape != b.shape:
+        return False
+
+    flat_a = a.flatten()
+    flat_b = b.flatten()
+
+    for i, (val_a, val_b) in enumerate(zip(flat_a, flat_b)):
+        # Handle NaN cases
+        if is_special_value(val_a) and is_special_value(val_b):
+            float_a = float(val_a)
+            float_b = float(val_b)
+            # Both NaN
+            if np.isnan(float_a) and np.isnan(float_b):
+                continue
+            # Both infinite with same sign
+            elif np.isinf(float_a) and np.isinf(float_b) and np.sign(float_a) == np.sign(float_b):
+                continue
+            else:
+                return False
+        elif is_special_value(val_a) or is_special_value(val_b):
+            return False
+        else:
+            try:
+                assert_quad_equal(val_a, val_b, rtol, atol)
+            except AssertionError:
+                return False
+
+    return True