Polar decomposition QDWH

GNUmakefile

@@ -508,6 +508,7 @@ ifneq ($(only_unit),1)
         src/gemmA.cc \
         src/gemmC.cc \
         src/geqrf.cc \
+        src/geqrf_qdwh_full.cc \
         src/gesv.cc \
         src/gesvMixed.cc \
         src/gesv_nopiv.cc \
@@ -534,6 +535,7 @@ ifneq ($(only_unit),1)
         src/hetrf.cc \
         src/hetrs.cc \
         src/norm.cc \
+        src/norm2est.cc \
         src/pbsv.cc \
         src/pbtrf.cc \
         src/pbtrs.cc \
@@ -543,6 +545,7 @@ ifneq ($(only_unit),1)
         src/potri.cc \
         src/potrs.cc \
         src/print.cc \
+        src/qdwh.cc \
         src/scale.cc \
         src/scale_row_col.cc \
         src/set.cc \
@@ -565,6 +568,7 @@ ifneq ($(only_unit),1)
         src/unmqr.cc \
         src/unmtr_hb2st.cc \
         src/unmtr_he2hb.cc \
+        src/unmqr_qdwh_full.cc \
         src/work/work_trmm.cc \
         src/work/work_trsm.cc \
         src/work/work_trsmA.cc \
@@ -672,6 +676,7 @@ ifneq ($(have_fortran),)
         test/pdlantr.f \
         test/pclantr.f \
         test/pzlantr.f \
+        test/test_qdwh.cc \
         # End. Add alphabetically, by base name after precision.
 endif
 endif

include/slate/slate.hh

@@ -1107,6 +1107,42 @@ void trcondest(
         blas::real_type<scalar_t> *rcond,
         Options const& opts = Options());
 
+//-----------------------------------------
+// Polar Decomposition
+
+//-----------------------------------------
+// qdwh()
+template <typename scalar_t>
+void qdwh(
+    Matrix<scalar_t>& A,
+    Matrix<scalar_t>& B,
+    int& itqr, int& itpo,
+    Options const& opts = Options());
+
+//-----------------------------------------
+// normest()
+template <typename matrix_type>
+blas::real_type<typename matrix_type::value_type>
+norm2est(
+    matrix_type &A,
+    Options const& opts = Options());
+
+//-----------------------------------------
+// geqrf_qdwh_full()
+template <typename scalar_t>
+void geqrf_qdwh_full(
+    Matrix<scalar_t>& A, TriangularFactors<scalar_t>& T,
+    Options const& opts = Options());
+
+//-----------------------------------------
+// unmqr_qdwh_full()
+template <typename scalar_t>
+void unmqr_qdwh_full(
+    Side side, Op op,
+    Matrix<scalar_t>& A, TriangularFactors<scalar_t>& T,
+    Matrix<scalar_t>& C,
+    Options const& opts = Options());
+
 } // namespace slate
 
 //-----------------------------------------