diff --git a/AMDiS/AMDISConfig.cmake.in b/AMDiS/AMDISConfig.cmake.in index aafcef469fd9996ee3910c37451a005f87d58188..0730cac200484df06fc77ace2f7dd87b188cd19d 100644 --- a/AMDiS/AMDISConfig.cmake.in +++ b/AMDiS/AMDISConfig.cmake.in @@ -108,6 +108,9 @@ if(AMDIS_NEED_BDDCML) set(AMDIS_BDDCML_LIB @BDDCML_LIB@) list(APPEND AMDIS_LIBRARIES ${AMDIS_BDDCML_LIB}) + + set(AMDIS_BDDCML_LINK_LIST @BDDCML_LINK_LIST@) + list(APPEND AMDIS_LIBRARIES ${AMDIS_BDDCML_LINK_LIST}) endif(AMDIS_NEED_BDDCML) diff --git a/AMDiS/CMakeLists.txt b/AMDiS/CMakeLists.txt index dc3019a72f377e8930706e9ca21f16e98263078c..2dfc618e85641bc8413e75121602e5cf6892e4f1 100644 --- a/AMDiS/CMakeLists.txt +++ b/AMDiS/CMakeLists.txt @@ -35,10 +35,9 @@ endif() SET(ENABLE_PARALLEL_DOMAIN "OFF" CACHE STRING "use parallel domain decomposition. please set to one of: PMTL, PETSC, OFF" ) option(USE_PETSC_DEV false) option(ENABLE_ZOLTAN false) -option(ENABLE_UMFPACK "use umfpack solver" false) +option(ENABLE_UMFPACK "Use of UMFPACK solver" false) option(ENABLE_BDDCML "Use of BDDCML library" false) - if(ENABLE_INTEL) Message("please set the icc manually") INCLUDE (CMakeForceCompiler) @@ -285,6 +284,8 @@ endif(ENABLE_UMFPACK) if(ENABLE_BDDCML) + SET(BDDCML_LINK_LIST "" CACHE STRING "Further libraries to be linked with BDDCML") + find_file(BDDCML_H bddcml_interface_c.h HINTS ENV CPATH DOC "Header bddcml_interface_c.h for the BDDCML library.") @@ -307,6 +308,7 @@ if(ENABLE_BDDCML) else() message(FATAL_ERROR "Could not find the BDDCML library") endif() + endif(ENABLE_BDDCML) diff --git a/AMDiS/src/parallel/BddcMlSolver.cc b/AMDiS/src/parallel/BddcMlSolver.cc index 260d7e5da6c12ee5b65e2724112a59097f8588e3..fba64a46f331fa55962b6eddde7d1b38f0c988ee 100644 --- a/AMDiS/src/parallel/BddcMlSolver.cc +++ b/AMDiS/src/parallel/BddcMlSolver.cc @@ -76,6 +76,9 @@ namespace AMDiS { int verboseLevel = 2; int numbase = 0; +// MSG("call to \"bddcml_init\" with the following arguments:\n"); +// MSG(" %d, [%d, %d], %d, %d, %d, %d, %d\n", nLevel, nSubdomains[0], nSubdomains[1], nLevel, nSubPerProc, c2f, verboseLevel, numbase); + bddcml_init(&nLevel, nSubdomains, &nLevel, &nSubPerProc, &c2f, &verboseLevel, &numbase); @@ -225,11 +228,56 @@ namespace AMDiS { // Number of non-zero entries in matrix int la = i_sparse.size(); - MSG("LOCAL LA = %d\n", la); + // MSG("LOCAL LA = %d\n", la); // Matrix is assembled int is_assembled_int = 0; + /* + string tmp=""; + + MSG("call to \"bddcml_upload_subdomain_data\" with the following arguments (each in one line):\n"); + MSG(" nelem = %d\n", nelem); + MSG(" nnod = %d\n", nnod); + MSG(" ndof = %d\n", ndof); + MSG(" ndim = %d\n", ndim); + MSG(" meshdim = %d\n", meshdim); + MSG(" isub = %d\n", isub); + MSG(" nelems = %d\n", nelems); + MSG(" nnods = %d\n", nnods); + MSG(" ndofs = %d\n", ndofs); + MSG(" inet = [%d, %d, %d, %d, %d, %d]\n", inet[0], inet[1], inet[2], inet[3], inet[4], inet[5]); + MSG(" linet = %d\n", linet); + MSG(" nnet = [%d, %d]\n", nnet[0], nnet[1]); + MSG(" lnnet = %d\n", nelems); + MSG(" nndf = [%d, %d, %d, %d]\n", nndf[0], nndf[1], nndf[2], nndf[3]); + MSG(" lnndf = %d\n", nnods); + MSG(" isngn = [%d, %d, %d, %d]\n", isngn[0], isngn[1], isngn[2], isngn[3]); + MSG(" lisngn = %d\n", nnods); + MSG(" isvgvn = [%d, %d, %d, %d]\n", isvgvn[0], isvgvn[1], isvgvn[2], isvgvn[3]); + MSG(" lisvgvn = %d\n", nnods); + MSG(" isegn = [%d, %d]\n", isegn[0], isegn[1]); + MSG(" lisegn = %d\n", nelems); + MSG(" xyz = [%f, %f, %f, %f, %f, %f, %f, %f]\n", xyz[0], xyz[1], xyz[2], xyz[3], xyz[4], xyz[5], xyz[6], xyz[7]); + MSG(" lxyz1 = %d\n", lxyz1); + MSG(" lxyz2 = %d\n", lxyz2); + MSG(" ifix = [%d, %d, %d, %d]\n", ifix[0], ifix[1], ifix[2], ifix[3]); + MSG(" lifix = %d\n", ndofs); + MSG(" fixv = [%f, %f, %f, %f]\n", fixv[0], fixv[1], fixv[2], fixv[3]); + MSG(" lfixv = %d\n", ndofs); + MSG(" rhs = [%f, %f, %f, %f]\n", rhs[0], rhs[1], rhs[2], rhs[3]); + MSG(" lrhs = %d\n", ndofs); + MSG(" is_rhs_complete = %d\n", is_rhs_complete); + MSG(" sol = [%f, %f, %f, %f]\n", sol[0], sol[1], sol[2], sol[3]); + MSG(" lsol = %d\n", ndofs); + MSG(" matrixtype = %d\n", matrixtype); + + MSG(" ispare = [%d, %d, %d, %d, %d, %d, %d, %d, %d, %d]\n", i_sparse[0], i_sparse[1], i_sparse[2], i_sparse[3], i_sparse[4], i_sparse[5], i_sparse[6], i_sparse[7], i_sparse[8], i_sparse[9]); + MSG(" jspare = [%d, %d, %d, %d, %d, %d, %d, %d, %d, %d]\n", j_sparse[0], j_sparse[1], j_sparse[2], j_sparse[3], j_sparse[4], j_sparse[5], j_sparse[6], j_sparse[7], j_sparse[8], j_sparse[9]); + MSG(" a_spare = [%f, %f, %f, %f, %f, %f, %f, %f, %f, %f]\n", a_sparse[0], a_sparse[1], a_sparse[2], a_sparse[3], a_sparse[4], a_sparse[5], a_sparse[6], a_sparse[7], a_sparse[8], a_sparse[9]); + MSG(" la = %d\n", la); + MSG(" is_assembled_int = %d\n", is_assembled_int); + */ bddcml_upload_subdomain_data(&nelem, &nnod, @@ -276,13 +324,18 @@ namespace AMDiS { int parallel_division_int = 1; int use_arithmetic_int = 1; int use_adaptive_int = 0; - MSG("BDDC POINT A\n"); + // MSG("call to \"bddcml_setup_preconditioner\" with the following arguments (each in one line):\n"); +// MSG(" %d\n", matrixtype); +// MSG(" %d\n", use_defaults_int); +// MSG(" %d\n", parallel_division_int); +// MSG(" %d\n", use_arithmetic_int); +// MSG(" %d\n", use_adaptive_int); + bddcml_setup_preconditioner(&matrixtype, &use_defaults_int, ¶llel_division_int, &use_arithmetic_int, &use_adaptive_int); - MSG("BDDC POINT B\n"); int method = 1; double tol = 1.e-6; @@ -301,8 +354,6 @@ namespace AMDiS { &converged_reason, &condition_number); - MSG("BDDC POINT C\n"); - MSG("BDDCML converged reason: %d within %d iterations \n", converged_reason, num_iter); @@ -354,8 +405,10 @@ namespace AMDiS { double val = value(*icursor); - i_sparse.push_back(rowIndex); - j_sparse.push_back(colIndex); +// i_sparse.push_back(rowIndex); +// j_sparse.push_back(colIndex); + i_sparse.push_back(*cursor * nComponents + ithRowComponent); + j_sparse.push_back(col(*icursor) * nComponents + ithColComponent); a_sparse.push_back(val); } } diff --git a/AMDiS/src/parallel/MeshDistributor.cc b/AMDiS/src/parallel/MeshDistributor.cc index 15bf0eef4e04a3d4398e35d53869488381060e5e..18a3388df420a165b77f5567421a9f476a5908f6 100644 --- a/AMDiS/src/parallel/MeshDistributor.cc +++ b/AMDiS/src/parallel/MeshDistributor.cc @@ -85,7 +85,8 @@ namespace AMDiS { repartitioningCounter(0), debugOutputDir(""), lastMeshChangeIndex(0), - createBoundaryDofFlag(0) + createBoundaryDofFlag(0), + sebastianMode(false) { FUNCNAME("MeshDistributor::ParalleDomainBase()"); @@ -238,10 +239,9 @@ namespace AMDiS { createInteriorBoundaryInfo(); -#if (DEBUG != 0) +#if (DEBUG != 0) ParallelDebug::printBoundaryInfo(*this); #endif - // === Remove neighbourhood relations due to periodic bounday conditions. === @@ -1576,8 +1576,18 @@ namespace AMDiS { bound.type = it->second; - AtomicBoundary& b = periodicBoundary.getNewAtomic(otherElementRank); - b = bound; + if (sebastianMode) { + if (bound.rankObj.elIndex == 3 && + bound.rankObj.ithObj == 1 || + bound.rankObj.elIndex == 78 && + bound.rankObj.ithObj == 2) { + AtomicBoundary& b = periodicBoundary.getNewAtomic(otherElementRank); + b = bound; + } + } else { + AtomicBoundary& b = periodicBoundary.getNewAtomic(otherElementRank); + b = bound; + } } } diff --git a/AMDiS/src/parallel/MeshDistributor.h b/AMDiS/src/parallel/MeshDistributor.h index cb73878554adbff6690ac015381afdb68b2a1568..0a028bda354ce414a9d1dd01f854800077e6fd20 100644 --- a/AMDiS/src/parallel/MeshDistributor.h +++ b/AMDiS/src/parallel/MeshDistributor.h @@ -711,6 +711,8 @@ namespace AMDiS { map<const FiniteElemSpace*, BoundaryDofInfo> boundaryDofInfo; public: + bool sebastianMode; + /// The boundary DOFs are sorted by subobject entities, i.e., first all /// face DOFs, edge DOFs and to the last vertex DOFs will be set to /// communication structure vectors, \ref sendDofs and \ref recvDofs. diff --git a/AMDiS/src/parallel/ParallelProblemStatBase.cc b/AMDiS/src/parallel/ParallelProblemStatBase.cc index 51d26281907e43e74188b075dc52a933d9311808..c7e4e8d7df4f008fee6b027f699dfeb5820b0243 100644 --- a/AMDiS/src/parallel/ParallelProblemStatBase.cc +++ b/AMDiS/src/parallel/ParallelProblemStatBase.cc @@ -31,6 +31,7 @@ namespace AMDiS { assembleMatrix, assembleVector); +#if (DEBUG != 0) double vm, rss; processMemUsage(vm, rss); vm /= 1024.0; @@ -42,6 +43,7 @@ namespace AMDiS { mpi::globalAdd(rss); MSG("Overall memory usage is VM = %.1f MB RSS = %.1f MB\n", vm, rss); +#endif } diff --git a/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc b/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc index 741970240a01d740990a2c378e0bbd334203fac1..5721a3e29baaa059a1153e0e7cb5644e6a1b34d9 100644 --- a/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc +++ b/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc @@ -25,7 +25,6 @@ namespace AMDiS { double wtime = MPI::Wtime(); vector<const FiniteElemSpace*> feSpaces = getFeSpaces(mat); - createGlobalDofMapping(feSpaces); int nRankRows = meshDistributor->getNumberRankDofs(feSpaces); int nOverallRows = meshDistributor->getNumberOverallDofs(feSpaces); @@ -43,6 +42,11 @@ namespace AMDiS { recvAllValues = 1; if (!d_nnz || recvAllValues != 0 || alwaysCreateNnzStructure) { + + // Global DOF mapping must only be recomputed, if the NNZ structure has + // changed (thus, also the mesh has changed). + createGlobalDofMapping(feSpaces); + if (d_nnz) { delete [] d_nnz; d_nnz = NULL; @@ -81,7 +85,7 @@ namespace AMDiS { for (int i = 0; i < nComponents; i++) for (int j = 0; j < nComponents; j++) if ((*mat)[i][j]) - setDofMatrix((*mat)[i][j], nComponents, i, j); + setDofMatrix((*mat)[i][j], i, j); #if (DEBUG != 0) MSG("Fill petsc matrix 2 needed %.5f seconds\n", MPI::Wtime() - wtime); @@ -122,7 +126,7 @@ namespace AMDiS { // === Transfer values from DOF vector to the PETSc vector. === for (int i = 0; i < vec->getSize(); i++) - setDofVector(petscRhsVec, vec->getDOFVector(i), vec->getSize(), i); + setDofVector(petscRhsVec, vec->getDOFVector(i), i); VecAssemblyBegin(petscRhsVec); VecAssemblyEnd(petscRhsVec); @@ -141,7 +145,7 @@ namespace AMDiS { VecSet(petscSolVec, 0.0); for (int i = 0; i < nComponents; i++) - setDofVector(petscSolVec, vec.getDOFVector(i), nComponents, i, true); + setDofVector(petscSolVec, vec.getDOFVector(i), i, true); VecAssemblyBegin(petscSolVec); VecAssemblyEnd(petscSolVec); @@ -193,8 +197,8 @@ namespace AMDiS { } - void PetscSolverGlobalMatrix::setDofMatrix(DOFMatrix* mat, int dispMult, - int dispAddRow, int dispAddCol) + void PetscSolverGlobalMatrix::setDofMatrix(DOFMatrix* mat, + int nRowMat, int nColMat) { FUNCNAME("PetscSolverGlobalMatrix::setDofMatrix()"); @@ -239,7 +243,7 @@ namespace AMDiS { // === Row DOF index is not periodic. === // Get PETSc's mat row index. - int rowIndex = dofToMatIndex.get(dispAddRow, globalRowDof); + int rowIndex = dofToMatIndex.get(nRowMat, globalRowDof); cols.clear(); values.clear(); @@ -253,7 +257,7 @@ namespace AMDiS { // Test if the current col dof is a periodic dof. bool periodicCol = perMap.isPeriodic(colFe, globalColDof); // Get PETSc's mat col index. - int colIndex = dofToMatIndex.get(dispAddCol, globalColDof); + int colIndex = dofToMatIndex.get(nColMat, globalColDof); // Ignore all zero entries, expect it is a diagonal entry. if (value(*icursor) == 0.0 && rowIndex != colIndex) @@ -302,7 +306,7 @@ namespace AMDiS { } for (unsigned int i = 0; i < newCols.size(); i++) { - cols.push_back(dofToMatIndex.get(dispAddCol, newCols[i])); + cols.push_back(dofToMatIndex.get(nColMat, newCols[i])); values.push_back(scaledValue); } } @@ -323,7 +327,6 @@ namespace AMDiS { // Traverse all column entries. for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); icursor != icend; ++icursor) { - // Global index of the current column index. int globalColDof = meshDistributor->mapDofToGlobal(colFe, col(*icursor)); @@ -392,8 +395,8 @@ namespace AMDiS { // === Translate the matrix entries to PETSc's matrix. for (unsigned int i = 0; i < entry.size(); i++) { - int rowIdx = dofToMatIndex.get(dispAddRow, entry[i].first); - int colIdx = dofToMatIndex.get(dispAddCol, entry[i].second); + int rowIdx = dofToMatIndex.get(nRowMat, entry[i].first); + int colIdx = dofToMatIndex.get(nColMat, entry[i].second); colsMap[rowIdx].push_back(colIdx); valsMap[rowIdx].push_back(scaledValue); @@ -419,8 +422,7 @@ namespace AMDiS { void PetscSolverGlobalMatrix::setDofVector(Vec& petscVec, DOFVector<double>* vec, - int dispMult, - int dispAdd, + int nRowVec, bool rankOnly) { FUNCNAME("PetscSolverGlobalMatrix::setDofVector()"); @@ -438,7 +440,7 @@ namespace AMDiS { DegreeOfFreedom globalRowDof = meshDistributor->mapDofToGlobal(feSpace, dofIt.getDOFIndex()); // Get PETSc's mat index of the row DOF. - int index = dofToMatIndex.get(dispAdd, globalRowDof); + int index = dofToMatIndex.get(nRowVec, globalRowDof); if (perMap.isPeriodic(feSpace, globalRowDof)) { std::set<int>& perAsc = perMap.getAssociations(feSpace, globalRowDof); @@ -448,7 +450,7 @@ namespace AMDiS { for (std::set<int>::iterator perIt = perAsc.begin(); perIt != perAsc.end(); ++perIt) { int mappedDof = perMap.map(feSpace, *perIt, globalRowDof); - int mappedIndex = dofToMatIndex.get(dispAdd, mappedDof); + int mappedIndex = dofToMatIndex.get(nRowVec, mappedDof); VecSetValues(petscVec, 1, &mappedIndex, &value, ADD_VALUES); } } else { diff --git a/AMDiS/src/parallel/PetscSolverGlobalMatrix.h b/AMDiS/src/parallel/PetscSolverGlobalMatrix.h index b9e827146a69770d225ef6112f77164035e584eb..e53add03bdef23aef666ec1eccf9d3e5d31f59d6 100644 --- a/AMDiS/src/parallel/PetscSolverGlobalMatrix.h +++ b/AMDiS/src/parallel/PetscSolverGlobalMatrix.h @@ -108,12 +108,11 @@ namespace AMDiS { void createPetscNnzStructure(Matrix<DOFMatrix*> *mat); /// Takes a DOF matrix and sends the values to the global PETSc matrix. - void setDofMatrix(DOFMatrix* mat, int dispMult = 1, - int dispAddRow = 0, int dispAddCol = 0); + void setDofMatrix(DOFMatrix* mat, int nRowMat = 0, int nColMat = 0); /// Takes a DOF vector and sends its values to a given PETSc vector. void setDofVector(Vec& petscVec, DOFVector<double>* vec, - int disMult = 1, int dispAdd = 0, bool rankOnly = false); + int nRowVec, bool rankOnly = false); void createGlobalDofMapping(vector<const FiniteElemSpace*> &feSpaces);