Skip to content

GitLab

Commit 4626656b authored by Thomas Witkowski's avatar Thomas Witkowski
Browse files

Okay, ab ins Wochenende.

parent 5e000286
Hide whitespace changes
Inline Side-by-side
Showing with 128 additions and 144 deletions
AMDiS/Reinit/src/HL_SignedDistTraverse.cc
View file @ 4626656b
......@@ -4,7 +4,7 @@
void HL_SignedDistTraverse::initializeBoundary()
{
FUNCNAME("HL_SignedDistTraverse::initializeBoundary()");
// ===== All non-boundary vertices are initialized with "infinity". =====
sD_DOF->set(inftyValue);
......
AMDiS/src/AdaptInfo.cc
View file @ 4626656b
......@@ -8,7 +8,7 @@ namespace AMDiS {
void AdaptInfo::setScalContents(int newSize)
{
int oldSize = scalContents.getSize();
int oldSize = static_cast<int>(scalContents.size());
if (newSize > oldSize) {
scalContents.resize(newSize);
......@@ -42,7 +42,7 @@ namespace AMDiS {
SerUtil::serialize(out, solverTolerance);
SerUtil::serialize(out, solverResidual);
int size = scalContents.getSize();
int size = static_cast<int>(scalContents.size());
SerUtil::serialize(out, size);
for (int i = 0; i < size; i++) {
SerUtil::serialize(out, scalContents[i]->est_sum);
......
AMDiS/src/AdaptInfo.h
View file @ 4626656b
......@@ -171,8 +171,9 @@ namespace AMDiS {
}
/// Destructor.
virtual ~AdaptInfo() {
for (int i = 0; i < scalContents.getSize(); i++)
virtual ~AdaptInfo()
{
for (unsigned int i = 0; i < scalContents.size(); i++)
delete scalContents[i];
}
......@@ -193,10 +194,10 @@ namespace AMDiS {
/// Returns whether space tolerance is reached.
virtual bool spaceToleranceReached()
{
int size = scalContents.getSize();
for (int i = 0; i < size; i++)
{
std::cout<<"est_sum:"<<scalContents[i]->est_sum<<" spaceTol:"<<scalContents[i]->spaceTolerance<<std::endl;
for (unsigned int i = 0; i < scalContents.size(); i++) {
std::cout << "est_sum:" <<scalContents[i]->est_sum
<< " spaceTol: " << scalContents[i]->spaceTolerance
<< std::endl;
if (!(scalContents[i]->est_sum < scalContents[i]->spaceTolerance))
return false;
}
......@@ -216,8 +217,7 @@ namespace AMDiS {
/// Returns whether time tolerance is reached.
virtual bool timeToleranceReached()
{
int size = scalContents.getSize();
for (int i = 0; i < size; i++)
for (unsigned int i = 0; i < scalContents.size(); i++)
if (!(scalContents[i]->est_t_sum < scalContents[i]->timeTolerance))
return false;
......@@ -236,8 +236,7 @@ namespace AMDiS {
/// Returns whether time error is under its lower bound.
virtual bool timeErrorLow()
{
int size = scalContents.getSize();
for (int i = 0; i < size; i++)
for (unsigned int i = 0; i < scalContents.size(); i++)
if (!(scalContents[i]->est_t_sum < scalContents[i]->timeErrLow))
return false;
......@@ -247,7 +246,7 @@ namespace AMDiS {
/// Print debug information about time error and its bound.
void printTimeErrorLowInfo()
{
for (int i = 0; i < scalContents.getSize(); i++)
for (unsigned int i = 0; i < scalContents.size(); i++)
std::cout << " Time error estimate = " << scalContents[i]->est_t_sum
<< " Time error bound = " << scalContents[i]->timeErrLow << "\n";
}
......@@ -387,6 +386,11 @@ namespace AMDiS {
/// Returns \ref est_sum.
inline double getEstSum(int index)
{
FUNCNAME("AdaptInfo::getEstSum()");
TEST_EXIT_DBG(static_cast<unsigned int>(index) < scalContents.size())
("Wrong index for adaptInfo!\n");
return scalContents[index]->est_sum;
}
......@@ -399,6 +403,11 @@ namespace AMDiS {
/// Returns \ref est_max.
inline double getEstMax(int index)
{
FUNCNAME("AdaptInfo::getEstSum()");
TEST_EXIT_DBG(static_cast<unsigned int>(index) < scalContents.size())
("Wrong index for adaptInfo!\n");
return scalContents[index]->est_max;
}
......@@ -586,7 +595,7 @@ namespace AMDiS {
inline int getSize()
{
return scalContents.getSize();
return scalContents.size();
}
inline void setSolverIterations(int it)
......@@ -728,7 +737,7 @@ namespace AMDiS {
double solverResidual;
/// Scalar adapt infos.
Vector<ScalContent*> scalContents;
std::vector<ScalContent*> scalContents;
/// Is true, if the adaptive procedure was deserialized from a file.
bool isDeserialized_;
......
AMDiS/src/Assembler.cc
View file @ 4626656b
......@@ -135,12 +135,12 @@ namespace AMDiS {
ElementMatrix &m =
smallElInfo->getSubElemCoordsMat(rowFESpace->getBasisFcts()->getDegree());
if (smallElInfo == colElInfo)
if (smallElInfo == colElInfo)
tmpMat = m * mat;
else
tmpMat = mat * trans(m);
mat = tmpMat;
tmpMat = mat * trans(m);
mat = tmpMat;
}
if (rememberElMat && &userMat != &elementMatrix)
......
AMDiS/src/DOFAdmin.cc
View file @ 4626656b
......@@ -246,14 +246,14 @@ namespace AMDiS {
ERROR("container not in list\n");
}
void DOFAdmin::compress(std::vector<DegreeOfFreedom> &new_dof)
{
FUNCNAME("DOFAdmin::compress()");
// nothing to do ?
if (size < 1) return;
if (usedCount < 1) return;
if (holeCount < 1) return;
if (size < 1 || usedCount < 1 || holeCount < 1)
return;
// vector to mark used dofs
for (int i = 0; i < size; i++)
......@@ -289,7 +289,7 @@ namespace AMDiS {
// get index of first changed dof
int first = last;
for (int i = 0; i<size; i++) {
if ((new_dof[i] < i) && (new_dof[i] >= 0)) {
if (new_dof[i] < i && new_dof[i] >= 0) {
first = i;
break;
}
......@@ -307,6 +307,7 @@ namespace AMDiS {
(*dc)->compressDOFContainer(n, new_dof);
}
void DOFAdmin::setNumberOfDOFs(int i, int v)
{
FUNCNAME("DOFAdmin::setNumberOfDOFs()");
......
AMDiS/src/DOFMatrix.cc
View file @ 4626656b
......@@ -345,6 +345,7 @@ namespace AMDiS {
addElementMatrix(elementMatrix, bound, mainElInfo, NULL);
}
void DOFMatrix::finishAssembling()
{
// call the operatos cleanup procedures
......@@ -353,6 +354,7 @@ namespace AMDiS {
(*it)->finishAssembling();
}
// Should work as before
Flag DOFMatrix::getAssembleFlag()
{
......@@ -364,6 +366,7 @@ namespace AMDiS {
return fillFlag;
}
void DOFMatrix::axpy(double a, const DOFMatrix& x, const DOFMatrix& y)
{
matrix+= a * x.matrix + y.matrix;
......@@ -375,6 +378,7 @@ namespace AMDiS {
matrix*= b;
}
void DOFMatrix::addOperator(Operator *op, double* factor, double* estFactor)
{
operators.push_back(op);
......@@ -382,6 +386,7 @@ namespace AMDiS {
operatorEstFactor.push_back(estFactor);
}
void DOFMatrix::serialize(std::ostream &out)
{
using namespace mtl;
......
AMDiS/src/DOFVector.cc
View file @ 4626656b
......@@ -829,6 +829,7 @@ namespace AMDiS {
addElementVector(factor, this->elementVector, bound, elInfo);
}
template<>
void DOFVectorBase<double>::assemble2(double factor,
ElInfo *mainElInfo, ElInfo *auxElInfo,
......
AMDiS/src/Element.cc
View file @ 4626656b
......@@ -260,7 +260,7 @@ namespace AMDiS {
void Element::newDOFFct2(const DOFAdmin* admin)
{
int i, j, k, n0, nd0;
int j, k, n0, nd0;
DegreeOfFreedom *ldof;
int vertices = mesh->getGeo(VERTEX);
int edges = mesh->getGeo(EDGE);
......@@ -270,7 +270,7 @@ namespace AMDiS {
if (nd) {
nd0 = admin->getNumberOfPreDOFs(VERTEX);
n0 = admin->getMesh()->getNode(VERTEX);
for (i = 0; i < vertices; i++) {
for (int i = 0; i < vertices; i++) {
CHANGE_DOFS_2(this);
}
}
......@@ -280,7 +280,7 @@ namespace AMDiS {
if (nd) {
nd0 = admin->getNumberOfPreDOFs(EDGE);
n0 = admin->getMesh()->getNode(EDGE);
for (i = 0; i < edges; i++) {
for (int i = 0; i < edges; i++) {
CHANGE_DOFS_2(this);
}
}
......@@ -291,7 +291,7 @@ namespace AMDiS {
if (nd) {
nd0 = admin->getNumberOfPreDOFs(FACE);
n0 = admin->getMesh()->getNode(FACE);
for (i = 0; i < faces; i++) {
for (int i = 0; i < faces; i++) {
CHANGE_DOFS_2(this);
}
}
......@@ -302,7 +302,7 @@ namespace AMDiS {
nd0 = admin->getNumberOfPreDOFs(CENTER);
n0 = admin->getMesh()->getNode(CENTER);
// only one center
i = 0;
int i = 0;
CHANGE_DOFS_2(this);
}
}
......
AMDiS/src/Marker.cc
View file @ 4626656b
......@@ -35,7 +35,7 @@ namespace AMDiS {
{
FUNCNAME("Marker::initMarking()");
TEST_EXIT(adaptInfo)("No AdaptInfo object!\n");
TEST_EXIT_DBG(adaptInfo)("No AdaptInfo object!\n");
elMarkRefine = 0;
elMarkCoarsen = 0;
......@@ -70,8 +70,12 @@ namespace AMDiS {
Flag Marker::markMesh(AdaptInfo *adaptInfo, Mesh *mesh)
{
FUNCNAME("Marker::markMesh()");
TEST_EXIT_DBG(mesh)("No mesh!\n");
initMarking(adaptInfo, mesh);
if (!adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row) &&
!adaptInfo->isRefinementAllowed(row == -1 ? 0 : row))
return 0;
......@@ -100,7 +104,7 @@ namespace AMDiS {
FUNCNAME("ESMarker::initMarking()");
Marker::initMarking(adaptInfo, mesh);
double ESThetaP = pow(ESTheta, p);
double ESThetaCP = pow(ESThetaC, p);
......
AMDiS/src/Mesh.cc
View file @ 4626656b
......@@ -63,7 +63,6 @@ namespace AMDiS {
// const Flag Mesh::USE_PARAMETRIC = 0X8000L ; // used in mg methods
DOFAdmin* Mesh::compressAdmin = NULL;
std::vector<DegreeOfFreedom> Mesh::dof_used;
const int Mesh::MAX_DOF = 100;
std::map<std::pair<DegreeOfFreedom, int>, DegreeOfFreedom*> Mesh::serializedDOFs;
......@@ -432,46 +431,39 @@ namespace AMDiS {
nNodeEl += 1;
}
void Mesh::dofCompress()
{
FUNCNAME("Mesh::dofCompress()");
Flag fill_flag;
for (int iadmin = 0; iadmin < static_cast<int>(admin.size()); iadmin++) {
compressAdmin = admin[iadmin];
for (unsigned int iadmin = 0; iadmin < admin.size(); iadmin++) {
DOFAdmin* compressAdmin = admin[iadmin];
TEST_EXIT_DBG(compressAdmin)("no admin[%d] in mesh\n", iadmin);
int size = compressAdmin->getSize();
if (size < 1)
continue;
if (compressAdmin->getUsedDOFs() < 1)
if (size < 1 ||
compressAdmin->getUsedDOFs() < 1 ||
compressAdmin->getHoleCount() < 1)
continue;
if (compressAdmin->getHoleCount() < 1)
continue;
newDOF.resize(size);
compressAdmin->compress(newDOF);
if (preserveCoarseDOFs)
fill_flag = Mesh::CALL_EVERY_EL_PREORDER | Mesh::FILL_NOTHING;
else
fill_flag = Mesh::CALL_LEAF_EL | Mesh::FILL_NOTHING;
Flag fill_flag = (preserveCoarseDOFs ?
Mesh::CALL_EVERY_EL_PREORDER | Mesh::FILL_NOTHING :
Mesh::CALL_LEAF_EL | Mesh::FILL_NOTHING);
TraverseStack stack;
ElInfo *elInfo = stack.traverseFirst(this, -1, fill_flag);
while (elInfo) {
newDOFFct1(elInfo);
elInfo->getElement()->newDOFFct1(compressAdmin);
elInfo = stack.traverseNext(elInfo);
}
elInfo = stack.traverseFirst(this, -1, fill_flag);
while (elInfo) {
newDOFFct2(elInfo);
elInfo->getElement()->newDOFFct2(compressAdmin);
elInfo = stack.traverseNext(elInfo);
}
......@@ -941,18 +933,6 @@ namespace AMDiS {
diam[i] = w;
}
int Mesh::newDOFFct1(ElInfo* ei)
{
ei->getElement()->newDOFFct1(compressAdmin);
return 0;
}
int Mesh::newDOFFct2(ElInfo* ei)
{
ei->getElement()->newDOFFct2(compressAdmin);
return 0;
}
void Mesh::serialize(std::ostream &out)
{
serializedDOFs.clear();
......
AMDiS/src/Mesh.h
View file @ 4626656b
......@@ -68,17 +68,6 @@ namespace AMDiS {
/// Assignment operator
Mesh& operator=(const Mesh&);
/** \name static methods used while mesh traversal
* \{
*/
/// Used while dof compress
static int newDOFFct1(ElInfo* e);
/// Used while dof compress
static int newDOFFct2(ElInfo* e);
/** \} */
/** \name getting methods
* \{
*/
......@@ -765,9 +754,6 @@ namespace AMDiS {
/// Needed during DOF compression (\ref DOFAdmin::compress).
std::vector<DegreeOfFreedom> newDOF;
/// Needed during DOF compression (\ref DOFAdmin::compress).
static DOFAdmin *compressAdmin;
/// Used by check functions
static std::vector<DegreeOfFreedom> dof_used;
......
AMDiS/src/ProblemVec.cc
View file @ 4626656b
......@@ -564,16 +564,21 @@ namespace AMDiS {
// here is reached already because of time adaption
allowFirstRefinement();
TEST_EXIT_DBG(static_cast<unsigned int>(nComponents == marker.size()))
("Wrong number of markers!\n");
Flag markFlag = 0;
for (int i = 0; i < nComponents; i++)
for (int i = 0; i < nComponents; i++) {
if (marker[i])
markFlag |= marker[i]->markMesh(adaptInfo, componentMeshes[i]);
else
WARNING("No marker for component %d\n", i);
WARNING("No marker for component %d\n", i);
}
return markFlag;
}
Flag ProblemVec::refineMesh(AdaptInfo *adaptInfo)
{
FUNCNAME("ProblemVec::refineMesh()");
......@@ -587,6 +592,7 @@ namespace AMDiS {
return refineFlag;
}
Flag ProblemVec::coarsenMesh(AdaptInfo *adaptInfo)
{
FUNCNAME("ProblemVec::coarsenMesh()");
......@@ -600,6 +606,7 @@ namespace AMDiS {
return coarsenFlag;
}
Flag ProblemVec::oneIteration(AdaptInfo *adaptInfo, Flag toDo)
{
FUNCNAME("ProblemVec::oneIteration()");
......@@ -620,6 +627,7 @@ namespace AMDiS {
return StandardProblemIteration::oneIteration(adaptInfo, toDo);
}
void ProblemVec::buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag,
bool asmMatrix, bool asmVector)
{
......@@ -633,13 +641,15 @@ namespace AMDiS {
// printOpenmpTraverseInfo(this, true);
// std::cout << "ElInfo = " << ElInfo::subElemMatrices.size() << std::endl;
for (unsigned int i = 0; i < meshes.size(); i++)
meshes[i]->dofCompress();
clock_t first = clock();
#ifdef _OPENMP
double wtime = omp_get_wtime();
#endif
for (int i = 0; i < static_cast<int>(meshes.size()); i++)
meshes[i]->dofCompress();
Flag assembleFlag =
flag |
......@@ -653,7 +663,8 @@ namespace AMDiS {
if (useGetBound)
assembleFlag |= Mesh::FILL_BOUND;
traverseInfo.updateStatus();
// Used to calculate the overall number of non zero entries.
......@@ -734,6 +745,7 @@ namespace AMDiS {
assembledMatrix[i][j] = true;
if (assembleMatrix)
matrix->finishInsertion();
......@@ -741,18 +753,18 @@ namespace AMDiS {
matrix->getBoundaryManager()->exitMatrix(matrix);
if (matrix)
nnz += matrix->getBaseMatrix().nnz();
nnz += matrix->getBaseMatrix().nnz();
}
// And now assemble boundary conditions on the vectors
assembleBoundaryConditions(rhs->getDOFVector(i),
solution->getDOFVector(i),
componentMeshes[i],
assembleFlag);
solution->getDOFVector(i),
componentMeshes[i],
assembleFlag);
}
solverMatrix.setMatrix(*systemMatrix);
createPrecon();
INFO(info, 8)("fillin of assembled matrix: %d\n", nnz);
......@@ -781,14 +793,16 @@ namespace AMDiS {
{
FUNCNAME("ProblemVec::dualAssemble()");
for (unsigned int i = 0; i < meshes.size(); i++)
meshes[i]->dofCompress();
clock_t first = clock();
#ifdef _OPENMP
double wtime = omp_get_wtime();
#endif
for (int i = 0; i < static_cast<int>(meshes.size()); i++)
meshes[i]->dofCompress();
Flag assembleFlag =
flag |
(*systemMatrix)[0][0]->getAssembleFlag() |
......@@ -801,7 +815,7 @@ namespace AMDiS {
if (useGetBound)
assembleFlag |= Mesh::FILL_BOUND;
traverseInfo.updateStatus();
// Used to calculate the overall number of non zero entries.
......@@ -887,11 +901,11 @@ namespace AMDiS {
bool cont = dualTraverse.traverseFirst(meshes[0], meshes[1], -1, -1,
assembleFlag, assembleFlag, dualElInfo);
while (cont) {
while (cont) {
bool newEl0 = (dualElInfo.rowElInfo->getElement()->getIndex() != oldElIndex0);
bool newEl1 = (dualElInfo.colElInfo->getElement()->getIndex() != oldElIndex1);
oldElIndex0 = dualElInfo.rowElInfo->getElement()->getIndex();
oldElIndex1 = dualElInfo.colElInfo->getElement()->getIndex();
oldElIndex1 = dualElInfo.colElInfo->getElement()->getIndex();
for (int i = 0; i < nComponents; i++) {
for (int j = 0; j < nComponents; j++) {
......@@ -901,64 +915,48 @@ namespace AMDiS {
continue;
if (traverseInfo.eqSpaces(i, j)) {
ElInfo *elInfo = NULL;
if (componentMeshes[i] == meshes[0] && newEl0)
elInfo = dualElInfo.rowElInfo;
if (componentMeshes[i] == meshes[1] && newEl1)
elInfo = dualElInfo.colElInfo;
if (elInfo != NULL) {
if (useGetBound)
basisFcts->getBound(elInfo, bound);
if (matrix) {
if (matrix)
matrix->assemble(1.0, elInfo, bound);
// Take the matrix boundary manager from the public matrix,
// but assemble the boundary conditions on the thread private matrix.
if (matrix->getBoundaryManager())
matrix->getBoundaryManager()->fillBoundaryConditions(elInfo, matrix);
}
if (i == j)
rhs->getDOFVector(i)->assemble(1.0, elInfo, bound);
}
if (traverseInfo.difAuxSpace(i) && i == j) {
ElInfo *mainElInfo, *auxElInfo;
if (traverseInfo.getRowFeSpace(i)->getMesh() == meshes[0]) {
mainElInfo = dualElInfo.rowElInfo;
auxElInfo = dualElInfo.colElInfo;
} else {
mainElInfo = dualElInfo.colElInfo;
auxElInfo = dualElInfo.rowElInfo;
}
if (useGetBound)
basisFcts->getBound(mainElInfo, bound);
rhs->getDOFVector(i)->assemble2(1.0, mainElInfo, auxElInfo,
dualElInfo.smallElInfo, dualElInfo.largeElInfo, bound);
ElInfo *mainElInfo, *auxElInfo;
if (traverseInfo.getRowFeSpace(i)->getMesh() == meshes[0]) {
mainElInfo = dualElInfo.rowElInfo;
auxElInfo = dualElInfo.colElInfo;
} else {
mainElInfo = dualElInfo.colElInfo;
auxElInfo = dualElInfo.rowElInfo;
}
if (traverseInfo.difAuxSpace(i, j) && matrix) {
ElInfo *mainElInfo, *auxElInfo;
if (traverseInfo.getRowFeSpace(i)->getMesh() == meshes[0]) {
mainElInfo = dualElInfo.rowElInfo;
auxElInfo = dualElInfo.colElInfo;
} else {