gepdf: Add a polar decomposition checking function

testing/testing_zcheck.c

@@ -1146,7 +1146,7 @@ int check_zlauum( run_arg_list_t *args, cham_uplo_t uplo, CHAM_desc_t *descA, CH
  *******************************************************************************
  */
 int check_zxxtrf( run_arg_list_t *args, cham_mtxtype_t mtxtype, cham_uplo_t uplo,
-                  CHAM_desc_t *descA, CHAM_desc_t *descLU )
+                   CHAM_desc_t *descA, CHAM_desc_t *descLU )
 {
     int info_local, info_global;
     int M = descA->m;
@@ -2147,4 +2147,66 @@ int check_zgepdf_qr( run_arg_list_t *args,
     return info_global;
 }
 
+/**
+ ********************************************************************************
+ *
+ * @ingroup testing
+ *
+ * @brief Checks if a Chameleon Polar Decomposition is correct.
+ *
+ * @Warning Check only the general case for now.
+ *
+ *******************************************************************************
+ *
+ * @param[in,out] descA
+ *          The descriptor of the original matrix, on exit the matrix is modified.
+ *
+ * @param[in] descU
+ *          The descriptor of the orthogonal polar factor of the decomposition.
+ *
+ * @param[in] descH
+ *          The descriptor of the symmetric/hermitian polar factor of the decomposition.
+ *
+ * @retval 0 successfull comparison
+ *
+ *******************************************************************************
+ */
+int check_zxxpd( run_arg_list_t *args,
+                 CHAM_desc_t *descA, CHAM_desc_t *descU, CHAM_desc_t *descH )
+{
+    int info_local, info_global;
+    double Anorm, Rnorm, result;
+    double eps = LAPACKE_dlamch_work('e');
+
+    /* Compute ||A|| */
+    Anorm = CHAMELEON_zlange_Tile( ChamFrobeniusNorm, descA );
+
+    /* R = A - U * H */
+    CHAMELEON_zgemm_Tile( ChamNoTrans, ChamNoTrans, 1., descU, descH, -1., descA );
+
+    /* Compute ||R|| */
+    Rnorm = CHAMELEON_zlange_Tile( ChamFrobeniusNorm, descA );
+
+    result = Rnorm / (Anorm * eps);
+    run_arg_add_double( args, "||A||",         Anorm );
+    run_arg_add_double( args, "||A-fact(A)||", Rnorm );
+
+    if ( isnan(result) || isinf(result) || (result > 60.0) ) {
+        info_local = 1;
+    }
+    else{
+        info_local = 0;
+    }
+
+    /* Broadcasts the result from the main processus */
+#if defined(CHAMELEON_USE_MPI)
+    MPI_Allreduce( &info_local, &info_global, 1, MPI_INT, MPI_MAX, MPI_COMM_WORLD );
+#else
+    info_global = info_local;
+#endif
+
+    (void)args;
+    return info_global;
+}
+
 #endif /* defined(CHAMELEON_SIMULATION) */

testing/testing_zcheck.h

@@ -62,12 +62,14 @@ static inline int check_zgeqrs        ( run_arg_list_t *args, cham_trans_t trans
 static inline int check_zgelqs        ( run_arg_list_t *args, cham_trans_t trans, CHAM_desc_t *descA, double Bnorm, CHAM_desc_t *descR ) { return 0; }
 static inline int check_zqc           ( run_arg_list_t *args, cham_side_t side, cham_trans_t trans, CHAM_desc_t *descC, CHAM_desc_t *descQ, CHAM_desc_t *descCC ) { return 0; }
 
-/* Matrix Generators */ 
+/* Matrix Generators */
 static inline int check_zrankk        ( run_arg_list_t *args, int K, CHAM_desc_t *descA ) { return 0; }
 
 /* Polar decomposition */
 static inline int check_zgepdf_qr     ( run_arg_list_t *args, CHAM_desc_t *descA1, CHAM_desc_t *descA2,
                                         CHAM_desc_t *descQ1, CHAM_desc_t *descQ2, CHAM_desc_t *descAF1 ) { return 0; }
+static inline int check_zxxpd         ( run_arg_list_t *args,
+               	                        CHAM_desc_t *descA, CHAM_desc_t *descU, CHAM_desc_t *descH ) { return 0; }
 
 #else /* !defined(CHAMELEON_SIMULATION) */
 
@@ -101,12 +103,14 @@ int check_zgeqrs        ( run_arg_list_t *args, cham_trans_t trans, CHAM_desc_t
 int check_zgelqs        ( run_arg_list_t *args, cham_trans_t trans, CHAM_desc_t *descA, double Bnorm, CHAM_desc_t *descR );
 int check_zqc           ( run_arg_list_t *args, cham_side_t side, cham_trans_t trans, CHAM_desc_t *descC, CHAM_desc_t *descQ, CHAM_desc_t *descCC );
 
-/* Matrix Generators */ 
+/* Matrix Generators */
 int check_zrankk        ( run_arg_list_t *args, int K, CHAM_desc_t *descA );
 
 /* Polar decomposition */
 int check_zgepdf_qr     ( run_arg_list_t *args, CHAM_desc_t *descA1, CHAM_desc_t *descA2,
                           CHAM_desc_t *descQ1, CHAM_desc_t *descQ2, CHAM_desc_t *descAF1 );
+int check_zxxpd         ( run_arg_list_t *args,
+                          CHAM_desc_t *descA, CHAM_desc_t *descU, CHAM_desc_t *descH );
 
 #endif /* defined(CHAMELEON_SIMULATION) */