Mentions légales du service

Skip to content
Snippets Groups Projects
Commit a0743aa5 authored by FELŠÖCI Marek's avatar FELŠÖCI Marek
Browse files

Force symmetric matrices and solvers

parent 6d9a081c
No related branches found
No related tags found
No related merge requests found
...@@ -33,12 +33,6 @@ extern int nbPts; ...@@ -33,12 +33,6 @@ extern int nbPts;
/*! \brief RHS for my tests */ /*! \brief RHS for my tests */
extern void *rhs; 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 */ /*! \brief Simple precision accuracy flag */
extern int simplePrec; extern int simplePrec;
......
...@@ -27,12 +27,6 @@ int printHelp() { ...@@ -27,12 +27,6 @@ int printHelp() {
" -height: Height of the cylinder. Default is 4.\n" " -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" " For FEM/BEM, the BEM is done on the outer surface of the cylinder (without the upper and lower disks).\n"
"\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" " Choice of the tested algorithm (no default):\n"
" -gemvhmat: gemv HMAT.\n" " -gemvhmat: gemv HMAT.\n"
" -solvehmat: solve HMAT.\n" " -solvehmat: solve HMAT.\n"
......
...@@ -409,7 +409,7 @@ HMAT_desc_t *HMAT_generate_matrix( hmat_interface_t *hi ) { ...@@ -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_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); hi->set_low_rank_epsilon(hmatrix, mpf_hmat_settings.epsilon);
...@@ -424,7 +424,7 @@ HMAT_desc_t *HMAT_generate_matrix( hmat_interface_t *hi ) { ...@@ -424,7 +424,7 @@ HMAT_desc_t *HMAT_generate_matrix( hmat_interface_t *hi ) {
ctx.progress = &mpf_hmat_settings.progress; ctx.progress = &mpf_hmat_settings.progress;
ctx.progress->user_data = NULL; ctx.progress->user_data = NULL;
ctx.factorization = hmat_factorization_none; ctx.factorization = hmat_factorization_none;
ctx.lower_symmetric = symMatSolver; ctx.lower_symmetric = 1;
ctx.compression = (hmat_compression_algorithm_t*) compression_ctx.output; ctx.compression = (hmat_compression_algorithm_t*) compression_ctx.output;
if (use_simple_assembly) { if (use_simple_assembly) {
ctx.simple_compute = interactionKernel; ctx.simple_compute = interactionKernel;
......
...@@ -15,14 +15,5 @@ int computeKernelBEM(double *coord1, double *coord2, int self, double complex *k ...@@ -15,14 +15,5 @@ int computeKernelBEM(double *coord1, double *coord2, int self, double complex *k
} else } else
*kernel = 1 / (4 * M_PI * r); *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; return 0;
} }
...@@ -15,8 +15,6 @@ int nbPtsPerLeaf = 50; ...@@ -15,8 +15,6 @@ int nbPtsPerLeaf = 50;
int nbAlgoRuns = 1; int nbAlgoRuns = 1;
char *sTYPE_ELEMENT_ARRAY_test_fembem[] = {NULL}; char *sTYPE_ELEMENT_ARRAY_test_fembem[] = {NULL};
int coupled = 0; int coupled = 0;
int symMatContent = 1;
int symMatSolver = 1;
int simplePrec = 0; int simplePrec = 0;
int complexALGO = 1; int complexALGO = 1;
ScalarType stype = DOUBLE_COMPLEX; ScalarType stype = DOUBLE_COMPLEX;
...@@ -143,27 +141,6 @@ int main(int argc, char **argv) { ...@@ -143,27 +141,6 @@ int main(int argc, char **argv) {
printf(" Setting lambda = %f (with %f points per wavelength)\n", lambda, ptsPerLambda) ; 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') --- */ /* --- Choix de l'arithmetique de calcul (default is '-z') --- */
if (MpfArgHasName(&argc, argv, 1, "-s") > 0) { if (MpfArgHasName(&argc, argv, 1, "-s") > 0) {
simplePrec = 1 ; simplePrec = 1 ;
......
...@@ -33,31 +33,13 @@ int testHMAT(double * relative_error) { ...@@ -33,31 +33,13 @@ int testHMAT(double * relative_error) {
// must choose one for him // must choose one for him
if(mpf_hmat_settings.factorization_type == hmat_factorization_none) { if(mpf_hmat_settings.factorization_type == hmat_factorization_none) {
if (use_hodlr) { if (use_hodlr) {
if(symMatSolver) mpf_hmat_settings.factorization_type = hmat_factorization_hodlrsym;
mpf_hmat_settings.factorization_type = hmat_factorization_hodlrsym;
else
mpf_hmat_settings.factorization_type = hmat_factorization_hodlr;
} else { } else {
if(symMatSolver) // use LDLT for real matrices
// use LDLT for real matrices mpf_hmat_settings.factorization_type = complexALGO ? hmat_factorization_llt : hmat_factorization_ldlt;
mpf_hmat_settings.factorization_type = complexALGO ? hmat_factorization_llt : hmat_factorization_ldlt;
else
mpf_hmat_settings.factorization_type = hmat_factorization_lu;
} }
} }
/* 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; HMAT_desc_t *hdesc;
hdesc = HMAT_generate_matrix( hi ); hdesc = HMAT_generate_matrix( hi );
hmat_matrix_t *hmatrix = hdesc->hmatrix; hmat_matrix_t *hmatrix = hdesc->hmatrix;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment