From 57914a8cf2cf9dfeab1c770e0780eb236b96d10d Mon Sep 17 00:00:00 2001 From: Pierre Jolivet Date: Tue, 9 Jul 2024 16:45:48 +0200 Subject: [PATCH 1/3] PETSc 3.21.3 --- 3rdparty/getall | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/3rdparty/getall b/3rdparty/getall index ce692a6ac..6fbe98ea1 100755 --- a/3rdparty/getall +++ b/3rdparty/getall @@ -171,10 +171,10 @@ download('parmmg','https://github.com/MmgTools/ParMmg/archive/f78f4e3e1019137787 'https://github.com/MmgTools', 'parmmg.zip', '9a76fef6be1ab66af3a813409a0388c9'); -download('PETSc','https://web.cels.anl.gov/projects/petsc/download/release-snapshots/petsc-3.21.2.tar.gz', +download('PETSc','https://web.cels.anl.gov/projects/petsc/download/release-snapshots/petsc-3.21.3.tar.gz', 'https://web.cels.anl.gov/projects/petsc/download/release-snapshots', - 'petsc-3.21.2.tar.gz', - '0bb68b348861a796ddfc2e94e3fcd02d'); + 'petsc-3.21.3.tar.gz', + 'c18cb846f97dffb9319d43e9017896b3'); download('htool','https://github.com/htool-ddm/htool/archive/946875d79d0036afb4dc2c0c13c165a607d830df.zip', 'https://github.com/htool-ddm/', From 8f3e54e5512746dcdb7e92a749a5bceb90ac0752 Mon Sep 17 00:00:00 2001 From: Pierre Jolivet Date: Sat, 13 Jul 2024 15:26:20 +0200 Subject: [PATCH 2/3] Use MATTRANSPOSEVIRTUAL instead of MATHERMITIANTRANSPOSEVIRTUAL --- plugin/mpi/PETSc-code.hpp | 83 ++++++++++++++++++++++++++++++++------- 1 file changed, 69 insertions(+), 14 deletions(-) diff --git a/plugin/mpi/PETSc-code.hpp b/plugin/mpi/PETSc-code.hpp index a18da50e4..36c6c7c50 100644 --- a/plugin/mpi/PETSc-code.hpp +++ b/plugin/mpi/PETSc-code.hpp @@ -918,12 +918,12 @@ namespace PETSc { MatAssemblyEnd(ptA->_petsc, MAT_FINAL_ASSEMBLY); } if (ptParent) { + Mat** mat; + PetscInt M, N; PetscBool assembled; + MatNestGetSubMats(ptParent->_petsc, &M, &N, &mat); MatAssembled(ptParent->_petsc, &assembled); if (!assembled) { - PetscInt M, N; - Mat** mat; - MatNestGetSubMats(ptParent->_petsc, &M, &N, &mat); PetscBool assemble = PETSC_TRUE; for (PetscInt i = 0; i < M && assemble; ++i) { for (PetscInt j = 0; j < N && assemble; ++j) { @@ -934,12 +934,51 @@ namespace PETSc { } } } - if (assemble) { - MatAssemblyBegin(ptParent->_petsc, MAT_FINAL_ASSEMBLY); - MatAssemblyEnd(ptParent->_petsc, MAT_FINAL_ASSEMBLY); + if (assemble) assembled = PETSC_TRUE; + } + if (PetscDefined(USE_COMPLEX)) { + for (PetscInt i = 0; i < M; ++i) { + for (PetscInt j = 0; j < N; ++j) { + if (mat[i][j]) { + MatType type; + PetscBool isType[2]; + MatGetType(mat[i][j], &type); + PetscStrcmp(type, MATHERMITIANTRANSPOSEVIRTUAL, isType); + if (!isType[0]) + PetscStrcmp(type, MATTRANSPOSEVIRTUAL, isType + 1); + else + isType[1] = PETSC_FALSE; + if (isType[0] || isType[1]) { + Mat C, D; + if (isType[0]) + MatHermitianTransposeGetMat(mat[i][j], &C); + else + MatTransposeGetMat(mat[i][j], &C); + if (C == ptA->_petsc) { + PetscReal norm; + MatDuplicate(C, MAT_COPY_VALUES, &D); + MatRealPart(D); + MatAXPY(D, -1.0, C, SAME_NONZERO_PATTERN); + MatNorm(D, NORM_INFINITY, &norm); + MatDestroy(&D); + if (norm < PETSC_MACHINE_EPSILON && isType[0]) { + MatCreateTranspose(C, &D); + MatNestSetSubMat(ptParent->_petsc, i, j, D); + MatDestroy(&D); + } + else if (norm > PETSC_MACHINE_EPSILON && isType[1]) { + MatCreateHermitianTranspose(C, &D); + MatNestSetSubMat(ptParent->_petsc, i, j, D); + MatDestroy(&D); + } + i = M, j = N; + } + } + } + } } } - else { + if (assembled) { MatAssemblyBegin(ptParent->_petsc, MAT_FINAL_ASSEMBLY); MatAssemblyEnd(ptParent->_petsc, MAT_FINAL_ASSEMBLY); } @@ -3265,7 +3304,7 @@ namespace PETSc { } void prepareConvert(Mat A, Mat* B) { MatType type; - PetscBool isType; + PetscBool isType[2]; Mat** mat; PetscInt M, N; MatNestGetSubMats(A, &M, &N, &mat); @@ -3275,13 +3314,23 @@ namespace PETSc { for (PetscInt j = 0; j < N; ++j) { if (mat[i][j]) { MatGetType(mat[i][j], &type); - PetscStrcmp(type, MATHERMITIANTRANSPOSEVIRTUAL, &isType); - if (isType) { + PetscStrcmp(type, MATHERMITIANTRANSPOSEVIRTUAL, isType); + if (PetscDefined(USE_COMPLEX) && !isType[0]) + PetscStrcmp(type, MATTRANSPOSEVIRTUAL, isType + 1); + else + isType[1] = PETSC_FALSE; + if (isType[0] || isType[1]) { b.emplace_back(std::make_pair(std::make_pair(i, j), Mat( ))); Mat D = mat[i][j]; Mat C; - MatHermitianTransposeGetMat(D, &b.back( ).second); - MatHermitianTranspose(b.back( ).second, MAT_INITIAL_MATRIX, &C); + if (isType[0]) { + MatHermitianTransposeGetMat(D, &b.back( ).second); + MatHermitianTranspose(b.back( ).second, MAT_INITIAL_MATRIX, &C); + } + else { + MatTransposeGetMat(D, &b.back( ).second); + MatTranspose(b.back( ).second, MAT_INITIAL_MATRIX, &C); + } PetscObjectReference((PetscObject)b.back().second); MatNestSetSubMat(A, i, j, C); MatDestroy(&C); @@ -4763,10 +4812,16 @@ namespace PETSc { PetscObjectTypeCompareAny((PetscObject)mat[i][j], &isType, MATMPIDENSE, MATSEQDENSE, ""); PetscInt n = 0, m = 0; if(!isType) { + Mat C = NULL; PetscStrcmp(type, MATHERMITIANTRANSPOSEVIRTUAL, &isType); - if(isType) { - Mat C; + if(isType) MatHermitianTransposeGetMat(mat[i][j], &C); + else if(PetscDefined(USE_COMPLEX)) { + PetscStrcmp(type, MATTRANSPOSEVIRTUAL, &isType); + if(isType) + MatTransposeGetMat(mat[i][j], &C); + } + if(C) { PetscObjectTypeCompareAny((PetscObject)C, &isType, MATMPIDENSE, MATSEQDENSE, ""); if(isType) MatGetSize(C, &m, &n); type = MATHERMITIANTRANSPOSEVIRTUAL; From 1ec138ed81b40c78a47a9b7c8a6f59a815ae8218 Mon Sep 17 00:00:00 2001 From: cmd8 <68232338+cmd8@users.noreply.github.com> Date: Wed, 17 Jul 2024 14:35:13 -0400 Subject: [PATCH 3/3] Adds support for `tgv = -3`, `tgv = -30` to eliminate only columns in matrices (#309) --- examples/tutorial/tgv-test.edp | 18 ++++++++++-------- src/femlib/HashMatrix.cpp | 14 +++++++------- src/fflib/problem.cpp | 8 ++++++++ 3 files changed, 25 insertions(+), 15 deletions(-) diff --git a/examples/tutorial/tgv-test.edp b/examples/tutorial/tgv-test.edp index ec7ad3ffa..72fddc795 100644 --- a/examples/tutorial/tgv-test.edp +++ b/examples/tutorial/tgv-test.edp @@ -1,15 +1,17 @@ load "msh3" meshL Th = segment(3); -real [int] vtgv = [ -20,-10,-2,-1,0, tgv]; +real [int] vtgv = [ -30,-20,-10,-3,-2,-1,0, tgv]; real [int] vsym = [ 0,1]; varf va(u,v)= int1d(Th)(u*v*4*9.)+on(1,2,u=0);// mul by 12 to have integer coef .. // label : 1 on first DoF O -// 2 on last DoF 2 -// tgv : -20 set to zero of row and colnum of Dof with label 1 or 2 -// tgv : -10 set to zero of row of Dof with label 1 or 2 -// tgv : -2 set to zero of row and colnum of Dof with label 1 or 2 and set one on diagonal term -// tgv : -1 set to zero of row of Dof with label 1 or 2 and set one on diagonal term -// tgv >=0 set to tgv value on diagonal term of Dof with label 1 or 2 +// 2 on last DoF 2 +// tgv : -30 set column of Dof with label 1 or 2 to zero +// tgv : -20 set row and column of Dof with label 1 or 2 to zero +// tgv : -10 set row of Dof with label 1 or 2 to zero +// tgv : -3 set column of Dof with label 1 or 2 to zero and set diagonal term to one +// tgv : -2 set row and column of Dof with label 1 or 2 to zero and set diagonal term to one +// tgv : -1 set row of Dof with label 1 or 2 to zero and set diagonal term to one +// tgv >=0 set diagonal term of Dof with label 1 or 2 to tgv value fespace Vh(Th,P1); int symj = 0; @@ -23,4 +25,4 @@ for [i,tgvi:vtgv] for [i,j,aij:A] F(i,j) = aij; cout << " sym= " << symj << " tgv = " << tgvi << " matrix= " << F << "\n\n\n"; -} \ No newline at end of file +} diff --git a/src/femlib/HashMatrix.cpp b/src/femlib/HashMatrix.cpp index e58b2e393..08b424dfc 100644 --- a/src/femlib/HashMatrix.cpp +++ b/src/femlib/HashMatrix.cpp @@ -1209,28 +1209,28 @@ void HashMatrix::SetBC(char *wbc,double ttgv) ntgv++; if(ttgv<0) { - + if( wbc[ii] ) { for (I k=p[ii];k 1.0e-10 && std::abs(ttgv+30.0) > 1.0e-10) + aij[k]=0;// put row to zero } } else operator()(ii,ii)=ttgv; } - if( std::abs(ttgv+2.0) < 1.0e-10 || std::abs(ttgv+20.0) < 1.0e-10) // remove also columm tgv == -2 ..... + if( std::abs(ttgv+2.0) < 1.0e-10 || std::abs(ttgv+20.0) < 1.0e-10 || std::abs(ttgv+3.0) < 1.0e-10 || std::abs(ttgv+30.0) < 1.0e-10) // remove also columm tgv == -2,-3 ..... { CSC(); for(I jj=0; jj< this->n; ++jj) if( wbc[jj] ) { for (I k=p[jj];k gg(buf); @@ -10061,6 +10063,8 @@ bool AssembleVarForm(Stack stack,const MMesh & Th,const FESpace1 & Uh,const FESp tgv = -10; else if (std::abs(tgv+2.0) < 1.0e-10) tgv = -20; + else if (std::abs(tgv+3.0) < 1.0e-10) + tgv = -30; } int Nbcomp=Vh.N; @@ -10218,6 +10222,8 @@ bool AssembleVarForm(Stack stack,const MMesh & Th,const FESpace1 & Uh,const FESp tgv = -10; else if (std::abs(tgv+2.0) < 1.0e-10) tgv = -20; + else if (std::abs(tgv+3.0) < 1.0e-10) + tgv = -30; } int Nbcomp=Vh.N; @@ -10383,6 +10389,8 @@ bool AssembleVarForm(Stack stack,const MMesh & Th,const FESpace1 & Uh,const FESp tgv = -10; else if (std::abs(tgv+2.0) < 1.0e-10) tgv = -20; + else if (std::abs(tgv+3.0) < 1.0e-10) + tgv = -30; } int Nbcomp=Vh.N;