ENH: initialised set up for qr implementation in umath_linalg.c.src

numpy/linalg/umath_linalg.c.src

@@ -3233,6 +3233,28 @@ static void
     set_fp_invalid_or_clear(error_occurred);
 }
 
+static void
+@TYPE@_qr(char **args, npy_intp const *dimensions, npy_intp const *steps,
+          void *NPY_UNUSED(func))
+{
+    GELSD_PARAMS_t params;
+    int error_occurred = get_fp_invalid_and_clear();
+    fortran_int n, m;
+    fortran_int excess;
+
+    INIT_OUTER_LOOP_7
+
+    m = (fortran_int)dimensions[0];
+    n = (fortran_int)dimensions[1];
+    excess = m - n;
+
+    if (init_@lapack_func@(&params, m, n, nrhs)) {
+        LINEARIZE_DATA_t a_in, q_out, r_out;
+
+    }
+
+}
+
 /**end repeat**/
 
 #pragma GCC diagnostic pop
@@ -3306,6 +3328,7 @@ GUFUNC_FUNC_ARRAY_REAL_COMPLEX(svd_N);
 GUFUNC_FUNC_ARRAY_REAL_COMPLEX(svd_S);
 GUFUNC_FUNC_ARRAY_REAL_COMPLEX(svd_A);
 GUFUNC_FUNC_ARRAY_REAL_COMPLEX(lstsq);
+GUFUNC_FUNC_ARRAY_REAL_COMPLEX(qr);
 GUFUNC_FUNC_ARRAY_EIG(eig);
 GUFUNC_FUNC_ARRAY_EIG(eigvals);
 
@@ -3379,6 +3402,13 @@ static char lstsq_types[] = {
     NPY_CDOUBLE, NPY_CDOUBLE, NPY_DOUBLE, NPY_CDOUBLE, NPY_DOUBLE, NPY_INT, NPY_DOUBLE,
 };
 
+static char qr_types[] = {
+    NPY_FLOAT,   NPY_FLOAT,   NPY_FLOAT,
+    NPY_DOUBLE,  NPY_DOUBLE, NPY_DOUBLE,
+    NPY_CFLOAT,  NPY_FLOAT,  NPY_FLOAT,
+    NPY_CDOUBLE, NPY_DOUBLE, NPY_DOUBLE,
+};
+
 typedef struct gufunc_descriptor_struct {
     char *name;
     char *signature;
@@ -3587,6 +3617,15 @@ GUFUNC_DESCRIPTOR_t gufunc_descriptors [] = {
         4, 3, 4,
         FUNC_ARRAY_NAME(lstsq),
         lstsq_types
+    },
+    {
+        "qr",
+        "(m,m),->(m,k),(k,m)",
+        "QR decomposition of last two dimensions and broadcast to the rest.\n"\
+        "The input matrix must be square.\n",
+        4, 1, 2,
+        FUNC_ARRAY_NAME(qr),
+        qr_types
     }
 };