From a0743aa5db3d6988585b3670ec433e3ccd5ee7f5 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?FEL=C5=A0=C3=96CI=20Marek?= <marek.felsoci@inria.fr>
Date: Sun, 5 Nov 2023 17:55:26 +0100
Subject: [PATCH] Force symmetric matrices and solvers

---
 include/main.h |  6 ------
 src/help.c     |  6 ------
 src/hmat.c     |  4 ++--
 src/kernel.c   |  9 ---------
 src/main.c     | 23 -----------------------
 src/testHMAT.c | 24 +++---------------------
 6 files changed, 5 insertions(+), 67 deletions(-)

diff --git a/include/main.h b/include/main.h
index 8507ad4..6dbaa7e 100644
--- a/include/main.h
+++ b/include/main.h
@@ -33,12 +33,6 @@ extern int nbPts;
 /*! \brief RHS for my tests */
 extern void *rhs;
 
-/*! \brief Symmetry flag for the content of the matrices (0=nosym, 1=sym, 2=sym pos def, 1 by default) */
-extern int symMatContent;
-
-/*! \brief Symmetry flag for the solver (true by default) */
-extern int symMatSolver;
-
 /*! \brief Simple precision accuracy flag */
 extern int simplePrec;
 
diff --git a/src/help.c b/src/help.c
index 40beca8..18c0fc7 100644
--- a/src/help.c
+++ b/src/help.c
@@ -27,12 +27,6 @@ int printHelp() {
          "     -height: Height of the cylinder. Default is 4.\n"
          "     For FEM/BEM, the BEM is done on the outer surface of the cylinder (without the upper and lower disks).\n"
          "\n"
-         "     Choice of the matrix symmetry (default is --sym):\n"
-         "     --sym: Symmetric matrices and solvers.\n"
-         "     --nosym: Non-Symmetric matrices and solvers (HMAT only).\n"
-         "     --symcont: Symmetric content in a non symmetric matrix and non symmetric solvers.\n"
-         "     --posdef: Symmetric positive definite matrices (requires --sym or --symcont).\n"
-         "\n"
          "     Choice of the tested algorithm (no default):\n"
          "     -gemvhmat: gemv HMAT.\n"
          "     -solvehmat: solve HMAT.\n"
diff --git a/src/hmat.c b/src/hmat.c
index e59cac8..21d72cd 100644
--- a/src/hmat.c
+++ b/src/hmat.c
@@ -409,7 +409,7 @@ HMAT_desc_t *HMAT_generate_matrix( hmat_interface_t *hi ) {
   }
 
   //  hmat_admissibility_t *admissibilityCondition = hmat_create_admissibility_standard(3.0);
-  hmat_matrix_t* hmatrix = hi->create_empty_hmatrix_admissibility(cluster_tree, cluster_tree, symMatSolver, admissibilityCondition);
+  hmat_matrix_t* hmatrix = hi->create_empty_hmatrix_admissibility(cluster_tree, cluster_tree, 1, admissibilityCondition);
 
   hi->set_low_rank_epsilon(hmatrix, mpf_hmat_settings.epsilon);
 
@@ -424,7 +424,7 @@ HMAT_desc_t *HMAT_generate_matrix( hmat_interface_t *hi ) {
   ctx.progress = &mpf_hmat_settings.progress;
   ctx.progress->user_data = NULL;
   ctx.factorization = hmat_factorization_none;
-  ctx.lower_symmetric = symMatSolver;
+  ctx.lower_symmetric = 1;
   ctx.compression = (hmat_compression_algorithm_t*) compression_ctx.output;
   if (use_simple_assembly) {
     ctx.simple_compute = interactionKernel;
diff --git a/src/kernel.c b/src/kernel.c
index 5415bd0..70ab6f2 100644
--- a/src/kernel.c
+++ b/src/kernel.c
@@ -15,14 +15,5 @@ int computeKernelBEM(double *coord1, double *coord2, int self, double complex *k
   } else
     *kernel = 1 / (4 * M_PI * r);
 
-  if (symMatContent==0) { // non symmetric case
-    // s = prod scal nu.r avec nu=(1,1,1)/sqrt(3)
-    double s = (v[0]+v[1]+v[2])/1.7320508075688772;
-    double fs=1+atan(s)/M_PI;
-    *kernel *= fs;
-  }
-  if (symMatContent==2 && self) { // symmetric positive definite case : add 'n.max(A)' times Identity
-    *kernel *= nbPts;
-  }
   return 0;
 }
diff --git a/src/main.c b/src/main.c
index 1855ad2..93e30eb 100644
--- a/src/main.c
+++ b/src/main.c
@@ -15,8 +15,6 @@ int nbPtsPerLeaf = 50;
 int nbAlgoRuns = 1;
 char *sTYPE_ELEMENT_ARRAY_test_fembem[] = {NULL};
 int coupled = 0;
-int symMatContent = 1;
-int symMatSolver = 1;
 int simplePrec = 0;
 int complexALGO = 1;
 ScalarType stype = DOUBLE_COMPLEX;
@@ -143,27 +141,6 @@ int main(int argc, char **argv) {
     printf("   Setting lambda = %f (with %f points per wavelength)\n", lambda, ptsPerLambda) ;
   }
 
-  /* --- Choice of the matrix symmetry (true by default) --- */
-  if (MpfArgHasName(&argc, argv, 1, "--sym") > 0) {
-    symMatSolver = 1 ;
-    symMatContent = 1 ;
-    printf("Testing: Symmetric matrices and solvers.\n") ;
-  }
-  if (MpfArgHasName(&argc, argv, 1, "--nosym") > 0) {
-    symMatSolver = 0 ;
-    symMatContent = 0 ;
-    printf("Testing: Non-Symmetric matrices and solvers.\n") ;
-  }
-  if (MpfArgHasName(&argc, argv, 1, "--symcont") > 0) {
-    symMatSolver = 0 ;
-    symMatContent = 1 ;
-    printf("Testing: Symmetric matrices and non symmetric solvers.\n") ;
-  }
-  if (symMatContent==1 && MpfArgHasName(&argc, argv, 1, "--posdef") > 0) {
-    symMatContent = 2 ;
-    printf("Testing: Positive definite matrices.\n") ;
-  }
-
   /* --- Choix de l'arithmetique de calcul (default is '-z') --- */
   if (MpfArgHasName(&argc, argv, 1, "-s") > 0) {
     simplePrec = 1 ;
diff --git a/src/testHMAT.c b/src/testHMAT.c
index 96ad02c..55e01f8 100644
--- a/src/testHMAT.c
+++ b/src/testHMAT.c
@@ -33,31 +33,13 @@ int testHMAT(double * relative_error) {
   // must choose one for him
   if(mpf_hmat_settings.factorization_type == hmat_factorization_none) {
     if (use_hodlr) {
-      if(symMatSolver)
-        mpf_hmat_settings.factorization_type = hmat_factorization_hodlrsym;
-      else
-        mpf_hmat_settings.factorization_type = hmat_factorization_hodlr;
+      mpf_hmat_settings.factorization_type = hmat_factorization_hodlrsym;
     } else {
-      if(symMatSolver)
-        // use LDLT for real matrices
-        mpf_hmat_settings.factorization_type = complexALGO ? hmat_factorization_llt : hmat_factorization_ldlt;
-      else
-        mpf_hmat_settings.factorization_type = hmat_factorization_lu;
+      // use LDLT for real matrices
+      mpf_hmat_settings.factorization_type = complexALGO ? hmat_factorization_llt : hmat_factorization_ldlt;
     }
   }
 
-  /* Check the coherency between the factorization chosen in MPF (with --hmat-lu/llt/ldlt/hodlr/hodlr-sym)
-     and the symmetry of the matrix (chosen with --sym/--nosym/--symcont) */
-  if (testedAlgo==_solveHMAT)
-    switch(mpf_hmat_settings.factorization_type){
-      case hmat_factorization_hodlrsym:
-      case hmat_factorization_ldlt:
-      case hmat_factorization_llt:ASSERTQ(symMatSolver); break;
-      // No 'hmat_factorization_hodlr' here, because it works on sym matrices
-      case hmat_factorization_lu:ASSERTQ(!symMatSolver); break;
-      default:break;
-    }
-
   HMAT_desc_t *hdesc;
   hdesc = HMAT_generate_matrix( hi );
   hmat_matrix_t *hmatrix = hdesc->hmatrix;
-- 
GitLab