diff --git a/AMDiS/src/Assembler.cc b/AMDiS/src/Assembler.cc index e94dc4b77f84f46fbc237201ec274c74dad8918d..b01d9f3c5da348b098828cad81a1648333c46df0 100644 --- a/AMDiS/src/Assembler.cc +++ b/AMDiS/src/Assembler.cc @@ -680,8 +680,6 @@ namespace AMDiS { void Pre0::calculateElementMatrix(const ElInfo *elInfo, ElementMatrix *mat) { - double *c = GET_MEMORY(double, 1); - if (firstCall) { q00 = Q00PsiPhi::provideQ00PsiPhi(owner->getRowFESpace()->getBasisFcts(), owner->getColFESpace()->getBasisFcts(), @@ -691,20 +689,22 @@ namespace AMDiS { firstCall = false; } - c[0] = 0.0; + // c[0] = 0.0; + double c = 0.0; int myRank = omp_get_thread_num(); + int size = static_cast<int>(terms[myRank].size()); - for (int i = 0; i < static_cast<int>( terms[myRank].size()); i++) { - (static_cast<ZeroOrderTerm*>((terms[myRank][i])))->getC(elInfo, 1, c); + for (int i = 0; i < size; i++) { + (static_cast<ZeroOrderTerm*>((terms[myRank][i])))->getC(elInfo, 1, &c); } - c[0] *= elInfo->getDet(); + c *= elInfo->getDet(); if (symmetric) { for (int i = 0; i < nRow; i++) { - (*mat)[i][i] += c[0] * q00->getValue(i,i); + (*mat)[i][i] += c * q00->getValue(i,i); for (int j = i + 1; j < nCol; j++) { - double val = c[0] * q00->getValue(i, j); + double val = c * q00->getValue(i, j); (*mat)[i][j] += val; (*mat)[j][i] += val; } @@ -712,16 +712,12 @@ namespace AMDiS { } else { for (int i = 0; i < nRow; i++) for (int j = 0; j < nCol; j++) - (*mat)[i][j] += c[0]*q00->getValue(i,j); + (*mat)[i][j] += c * q00->getValue(i, j); } - - FREE_MEMORY(c, double, 1); } void Pre0::calculateElementVector(const ElInfo *elInfo, ElementVector *vec) { - double *c = GET_MEMORY(double, 1); - if (firstCall) { q00 = Q00PsiPhi::provideQ00PsiPhi(owner->getRowFESpace()->getBasisFcts(), owner->getColFESpace()->getBasisFcts(), @@ -734,17 +730,15 @@ namespace AMDiS { ::std::vector<OperatorTerm*>::iterator termIt; int myRank = omp_get_thread_num(); - c[0] = 0.0; + double c = 0.0; for (termIt = terms[myRank].begin(); termIt != terms[myRank].end(); ++termIt) { - (static_cast<ZeroOrderTerm*>( *termIt))->getC(elInfo, 1, c); + (static_cast<ZeroOrderTerm*>( *termIt))->getC(elInfo, 1, &c); } - c[0] *= elInfo->getDet(); + c *= elInfo->getDet(); for (int i = 0; i < nRow; i++) - (*vec)[i] += c[0] * q0->getValue(i); - - FREE_MEMORY(c, double, 1); + (*vec)[i] += c * q0->getValue(i); } Stand10::Stand10(Operator *op, Assembler *assembler, Quadrature *quad) diff --git a/AMDiS/src/DOFMatrix.cc b/AMDiS/src/DOFMatrix.cc index 413757fdb1ebdb6154083ac84b8dd65d92bd5060..4198906a9b44fdd1a32e18ce65818aa35b5afdfd 100644 --- a/AMDiS/src/DOFMatrix.cc +++ b/AMDiS/src/DOFMatrix.cc @@ -35,7 +35,7 @@ namespace AMDiS { { TEST_EXIT(rowFESpace)("no rowFESpace\n"); - if(!colFESpace) { + if (!colFESpace) { colFESpace = rowFESpace; } @@ -57,7 +57,7 @@ namespace AMDiS { { FUNCNAME("DOFMatrix::~DOFMatrix()"); - if(rowFESpace && rowFESpace->getAdmin()) { + if (rowFESpace && rowFESpace->getAdmin()) { (const_cast<DOFAdmin*>(rowFESpace->getAdmin()))->removeDOFIndexed(this); } } diff --git a/AMDiS/src/ElInfo1d.cc b/AMDiS/src/ElInfo1d.cc index c57751e9eb79cb1ec49fb537141a24047e431a75..7eceb7ece8b38e615bbb0b5a1428aeff86650a5c 100644 --- a/AMDiS/src/ElInfo1d.cc +++ b/AMDiS/src/ElInfo1d.cc @@ -33,11 +33,7 @@ namespace AMDiS { for (int i = 0; i < vertices; i++) { coord_[i] = mel->coord[i]; } - } - - // if(fillFlag_.isSet(Mesh::FILL_DET) || fillFlag_.isSet(Mesh::FILL_GRD_LAMBDA)) { - // det = elGrdLambda(*Lambda); - // } + } if (fillFlag_.isSet(Mesh::FILL_NEIGH) || fillFlag_.isSet(Mesh::FILL_OPP_COORDS)) { WorldVector<double> oppC; @@ -221,10 +217,6 @@ namespace AMDiS { } } - // if(fillFlag_.isSet(Mesh::FILL_DET) || fillFlag_.isSet(Mesh::FILL_GRD_LAMBDA)) { - // det = calcGrdLambda(*Lambda); - // } - if (fillFlag_.isSet(Mesh::FILL_NEIGH) || fillFlag_.isSet(Mesh::FILL_OPP_COORDS)) { WorldVector<double> oppC; diff --git a/AMDiS/src/ElInfo3d.cc b/AMDiS/src/ElInfo3d.cc index a8d89eb51d4a9c31786601999f22d9f596429811..5d71ba52d98feb10128a293a8372c6c927eb9971 100644 --- a/AMDiS/src/ElInfo3d.cc +++ b/AMDiS/src/ElInfo3d.cc @@ -348,7 +348,6 @@ namespace AMDiS { { FUNCNAME("ElInfo3d::fillElInfo()"); - int el_type_local = 0; /* el_type in {0,1,2} */ int ochild = 0; /* index of other child = 1-ichild */ int *cv = NULL; /* cv = child_vertex[el_type][ichild] */ const int (*cvg)[4] = NULL; /* cvg = child_vertex[el_type] */ @@ -370,12 +369,12 @@ namespace AMDiS { fillFlag_ = fillFlag__local; parent_ = el_old; level_ = elinfo_old->level_ + 1; - el_type = (( dynamic_cast<ElInfo3d*>(const_cast<ElInfo*>( elinfo_old)))->el_type + 1) % 3; + int el_type_local = ( dynamic_cast<ElInfo3d*>(const_cast<ElInfo*>( elinfo_old)))->getType(); + el_type = (el_type_local + 1) % 3; TEST_EXIT_DBG(element_)("missing child %d?\n", ichild); if (fillFlag__local.isAnySet()) { - el_type_local = ( dynamic_cast<ElInfo3d*>(const_cast<ElInfo*>( elinfo_old)))->getType(); cvg = Tetrahedron::childVertex[el_type_local]; cv = const_cast<int*>( cvg[ichild]); ochild = 1 - ichild; diff --git a/AMDiS/src/ElInfo3d.h b/AMDiS/src/ElInfo3d.h index 6546270ca2084383bf4f8e8e031182028d785aaf..eee554a7f2b3684f8866355a6797eb1d317dbf2d 100644 --- a/AMDiS/src/ElInfo3d.h +++ b/AMDiS/src/ElInfo3d.h @@ -81,13 +81,6 @@ namespace AMDiS { */ double getNormal(int side, WorldVector<double> &normal) const; - /** \brief - * 3-dimensional realisation of ElInfo's getElementNormal method. - */ - //double getElementNormal( WorldVector<double> &normal) const; - - - /** \brief * update ElInfo after refinement (of some neighbours). Only in 3d! */ diff --git a/AMDiS/src/FileWriter.cc b/AMDiS/src/FileWriter.cc index 1c9b9147ce770e41c3b6677efbd1088475f2f112..592292385f9430e75e0322dba14adcfbb0930459 100644 --- a/AMDiS/src/FileWriter.cc +++ b/AMDiS/src/FileWriter.cc @@ -10,6 +10,7 @@ #include "Flag.h" #include "ElInfo.h" #include "Mesh.h" +#include "OpenMP.h" namespace AMDiS { @@ -150,6 +151,10 @@ namespace AMDiS { GET_PARAMETER(0, name + "->index decimals", "%d", &indexDecimals); GET_PARAMETER(0, name + "->write every i-th timestep", "%d", &tsModulo); GET_PARAMETER(0, name + "->delay", "%d", &delayWriting_); + + TEST_EXIT(!delayWriting_ || amdisHaveOpenMP) + ("Delayed writing only possible with OpenMP support!\n"); + } void FileWriter::writeFiles(AdaptInfo *adaptInfo, diff --git a/AMDiS/src/Global.cc b/AMDiS/src/Global.cc index 84a4b3ff57b63caf385d2bd342107a87b90bd72a..1b80f9f491ca333c7a3d8e8cf608386dccf69d28 100644 --- a/AMDiS/src/Global.cc +++ b/AMDiS/src/Global.cc @@ -13,14 +13,18 @@ namespace AMDiS { const char *funcName = NULL; const char *Msg::oldFuncName = NULL; - ::std::ostream* Msg::out=NULL; - ::std::ostream* Msg::error=NULL; + ::std::ostream* Msg::out = NULL; + ::std::ostream* Msg::error = NULL; int Global::dimOfWorld = 0; + ::std::vector< ::std::vector< int > > Global::geoIndexTable; int Msg::msgInfo = 10; bool Msg::msgWait = true; Element *Global::referenceElement[4] = - { NULL, NEW Line(NULL), NEW Triangle(NULL), NEW Tetrahedron(NULL) }; + { NULL, + NEW Line(NULL), + NEW Triangle(NULL), + NEW Tetrahedron(NULL) }; void Msg::wait(bool w) { @@ -55,99 +59,87 @@ namespace AMDiS { void Msg::open_file(const char *filename, OPENMODE type) { - FUNCNAME("Msg::open_file"); + FUNCNAME("Msg::open_file()"); ::std::ofstream *fp; - if (filename && (fp = new ::std::ofstream(filename, type))) - { - if (out && *out != ::std::cout && *out != ::std::cerr) - { - dynamic_cast< ::std::ofstream*>(out)->close(); - delete out; - } - - out = fp; - } - else - { - if (filename) - ERROR("can not open %s;\n", filename); - else - ERROR("no filename specified;\n"); - ERROR("use previous stream for messages furthermore\n"); + if (filename && (fp = new ::std::ofstream(filename, type))) { + if (out && *out != ::std::cout && *out != ::std::cerr) { + dynamic_cast< ::std::ofstream*>(out)->close(); + delete out; } + + out = fp; + } else { + if (filename) + ERROR("can not open %s;\n", filename); + else + ERROR("no filename specified;\n"); + + ERROR("use previous stream for messages furthermore\n"); + } return; } void Msg::change_error_out(::std::ofstream *fp) { - FUNCNAME("Msg::change_error_out"); - if (fp) - { - if (error && *error != ::std::cout && *error != ::std::cerr) - { - dynamic_cast< ::std::ofstream*>(error)->close(); - delete error; - } - - error = fp; - } - else - { - ERROR("file pointer is pointer to nil;\n"); - ERROR("use previous stream for errors furthermore\n"); + FUNCNAME("Msg::change_error_out()"); + + if (fp) { + if (error && *error != ::std::cout && *error != ::std::cerr) { + dynamic_cast< ::std::ofstream*>(error)->close(); + delete error; } + + error = fp; + } else { + ERROR("file pointer is pointer to nil;\n"); + ERROR("use previous stream for errors furthermore\n"); + } return; } void Msg::open_error_file(const char *filename, OPENMODE type) { - FUNCNAME("Msg::open_error_file"); + FUNCNAME("Msg::open_error_file()"); ::std::ofstream *fp; - if (filename && (fp = new ::std::ofstream(filename, type))) - { - if (error && *error != ::std::cout && *error != ::std::cerr) - { - dynamic_cast< ::std::ofstream*>(error)->close(); - delete error; - } - - error = fp; - } - else - { - if (filename) - ERROR("can not open %s;\n", filename); - else - ERROR("no filename specified;\n"); - ERROR("use previous stream for errors furthermore\n"); + if (filename && (fp = new ::std::ofstream(filename, type))) { + if (error && *error != ::std::cout && *error != ::std::cerr) { + dynamic_cast< ::std::ofstream*>(error)->close(); + delete error; } + + error = fp; + } else { + if (filename) + ERROR("can not open %s;\n", filename); + else + ERROR("no filename specified;\n"); + ERROR("use previous stream for errors furthermore\n"); + } + return; } const char *generate_filename(const char * path, const char * fn, int ntime) { - static char name[256]; - char *cp; + static char name[256]; + char *cp; - if (path == NULL || path[0] == '\0') - { - sprintf(name, "./%s", fn); - } - else - { - const char* ccp = path; - while (*ccp) - ccp++; - ccp--; - if (*ccp == '/') - sprintf(name, "%s%s", path, fn); - else - sprintf(name, "%s/%s", path, fn); - } + if (path == NULL || path[0] == '\0') { + sprintf(name, "./%s", fn); + } else { + const char* ccp = path; + while (*ccp) + ccp++; + ccp--; + if (*ccp == '/') + sprintf(name, "%s%s", path, fn); + else + sprintf(name, "%s/%s", path, fn); + } cp = name; while (*cp) cp++; @@ -159,7 +151,8 @@ namespace AMDiS { void Msg::print_funcname(const char *funcName) { - if (!out) out = &::std::cout; + if (!out) + out = &::std::cout; if (funcName && oldFuncName != funcName) { (*out)<< funcName << ":" << ::std::endl; @@ -176,29 +169,32 @@ namespace AMDiS { { static int old_line = -1; - if (!error) error = &::std::cerr; + if (!error) + error = &::std::cerr; if (funcName && oldFuncName != funcName) { (*error)<<funcName<< ": "; } else if (!funcName) { - if (line-old_line > 5) (*error)<< "*unknown function*"; + if (line-old_line > 5) + (*error) << "*unknown function*"; } if (oldFuncName != funcName) { - (*error)<<"ERROR in "<<file<<", line "<<line<<"\n"; + (*error) << "ERROR in " << file << ", line " << line << "\n"; oldFuncName = funcName; } else if (line - old_line > 5) - (*error)<< "ERROR in "<<file<<", line "<<line<<"\n"; + (*error) << "ERROR in " << file << ", line " << line << "\n"; old_line = line; } void Msg::print_error_exit(const char *format, ...) { - va_list arg; - char buff[255]; + va_list arg; + char buff[255]; - if (!error) error = &::std::cerr; + if (!error) + error = &::std::cerr; va_start(arg, format); vsprintf(buff, format, arg); @@ -210,11 +206,12 @@ namespace AMDiS { void Msg::print_error(const char *format, ...) { - va_list arg; - char buff[255]; + va_list arg; + char buff[255]; - if (!error) error = &::std::cerr; + if (!error) + error = &::std::cerr; va_start(arg, format); vsprintf(buff, format, arg); @@ -228,31 +225,33 @@ namespace AMDiS { const char *file, int line) { - static int old_line = -1; + static int old_line = -1; - if (!out) out = &::std::cout; + if (!out) + out = &::std::cout; if (funcName && oldFuncName != funcName) { - (*out)<<funcName<<": "; + (*out) << funcName << ": "; } else if (!funcName){ - (*out)<< "*unknown function*"; + (*out) << "*unknown function*"; } if (oldFuncName != funcName) { - (*out)<<"WARNING in "<<file<<", line "<<line<<"\n"; + (*out) << "WARNING in " << file << ", line " << line << "\n"; oldFuncName = funcName; } else if (line - old_line > 5) - (*out)<<"WARNING in "<<file<<", line "<<line<<"\n"; + (*out) << "WARNING in " << file << ", line " << line << "\n"; old_line = line; } void Msg::print_warn(const char *format, ...) { - va_list arg; - char buff[255]; + va_list arg; + char buff[255]; - if (!out) out = &::std::cout; + if (!out) + out = &::std::cout; va_start(arg, format); vsprintf(buff, format, arg); @@ -265,15 +264,15 @@ namespace AMDiS { void Msg::print(const char *format, ...) { - va_list arg; - char buff[255]; + va_list arg; + char buff[255]; if (!out) out = &::std::cout; va_start(arg, format); vsprintf(buff, format, arg); - (*out)<<buff; + (*out) << buff; va_end(arg); } @@ -283,54 +282,53 @@ namespace AMDiS { // get dimension TEST_EXIT(Parameters::initialized())("Parameters not initialized!\n"); - Parameters::getGlobalParameter(0,"dimension of world","%d",&d); - TEST_EXIT(d > 0)("can not initialize dimension\n"); + Parameters::getGlobalParameter(0, "dimension of world","%d",&d); + TEST_EXIT(d > 0)("Cannot initialize dimension!\n"); TEST_EXIT((d == 1) || (d == 2) || (d == 3))("Invalid world dimension %d!\n",d); // set dimension dimOfWorld = d; - // set msgWait - Parameters::getGlobalParameter(0,"WAIT","%d",&d); - Msg::setMsgWait(!(d==0)); - }; - - int Global::getGeo(GeoIndex p, int dim) - { - FUNCNAME("Global::getGeo()"); - - initTest(); - TEST_EXIT_DBG((p >= MINPART) && (p <= MAXPART)) - ("Calling for invalid geometry value %d\n",p); - TEST_EXIT_DBG((dim >= 0) && (dim < 4)) - ("invalid dim: %d\n", dim); - - if (p == WORLD) - return dimOfWorld; - - if (dim == 0) { - switch(p) { - case PARTS: - return 1; - case VERTEX: - return 1; - case EDGE: - return 0; - case FACE: - return 0; - default: - ERROR_EXIT("dim = 0\n"); + // prepare geoIndex-Table + int geoTableSize = abs(static_cast<int>(MINPART)) + MAXPART + 1; + geoIndexTable.resize(4); + for (int i = 0; i < 4; i++) { + geoIndexTable[i].resize(geoTableSize); + for (int j = 0; j < geoTableSize; j++) { + geoIndexTable[i][j] = 0; } } - return referenceElement[dim]->getGeo(p); - } + geoIndexTable[0][PARTS - MINPART] = 1; + geoIndexTable[0][VERTEX - MINPART] = 1; + geoIndexTable[0][EDGE - MINPART] = 0; + geoIndexTable[0][FACE - MINPART] = 0; + geoIndexTable[0][WORLD - MINPART] = dimOfWorld; + + for (int i = 1; i < 4; i++) { + geoIndexTable[i][CENTER - MINPART] = referenceElement[i]->getGeo(CENTER); + geoIndexTable[i][VERTEX - MINPART] = referenceElement[i]->getGeo(VERTEX); + geoIndexTable[i][EDGE - MINPART] = referenceElement[i]->getGeo(EDGE); + geoIndexTable[i][FACE - MINPART] = referenceElement[i]->getGeo(FACE); + geoIndexTable[i][DIMEN - MINPART] = referenceElement[i]->getGeo(DIMEN); + geoIndexTable[i][PARTS - MINPART] = referenceElement[i]->getGeo(PARTS); + geoIndexTable[i][NEIGH - MINPART] = referenceElement[i]->getGeo(NEIGH); + geoIndexTable[i][WORLD - MINPART] = dimOfWorld; + geoIndexTable[i][BOUNDARY - MINPART] = referenceElement[i]->getGeo(BOUNDARY); + geoIndexTable[i][PROJECTION - MINPART] = referenceElement[i]->getGeo(PROJECTION); + } + // set msgWait + Parameters::getGlobalParameter(0, "WAIT", "%d", &d); + Msg::setMsgWait(!(d == 0)); + }; int fac(int i) { - if(i<=1) return 1; - else return i*fac(i-1); + if (i <= 1) + return 1; + else + return i * fac(i - 1); } ::std::string memSizeStr(int size) diff --git a/AMDiS/src/Global.h b/AMDiS/src/Global.h index 2dc0a706d75436842aab5477d960c60181fe7ee7..9cf7fae19035f85094d9087e779e2732cd968790 100644 --- a/AMDiS/src/Global.h +++ b/AMDiS/src/Global.h @@ -453,8 +453,6 @@ namespace AMDiS { * returns geometrical information. Currently this is only dimOfWorld. */ static inline int getGeo(GeoIndex p) { - FUNCNAME("Global::getGeo()"); - initTest(); if (WORLD == p) return dimOfWorld; @@ -467,7 +465,17 @@ namespace AMDiS { * returns geometrical information about elements of the dimension dim. * getGeo(VERTEX, 3) returns 4 because a Tetrahedron has 4 vertices. */ - static int getGeo(GeoIndex p, int dim); + static inline int getGeo(GeoIndex p, int dim) { + initTest(); + TEST_EXIT_DBG((p >= MINPART) && (p <= MAXPART)) + ("Calling for invalid geometry value %d\n",p); + TEST_EXIT_DBG((dim >= 0) && (dim < 4)) + ("invalid dim: %d\n", dim); + TEST_EXIT_DBG((dim != 0) || (p == PARTS || p == VERTEX || p == EDGE || p == FACE)) + ("dim = 0\n"); + + return geoIndexTable[dim][p - MINPART]; + } private: /** \brief @@ -503,6 +511,8 @@ namespace AMDiS { * Tetrahedron wich is 4 => no switch statement is necessary. */ static Element *referenceElement[4]; + + static ::std::vector< ::std::vector< int > > geoIndexTable; }; #define COMMA , diff --git a/AMDiS/src/Lagrange.cc b/AMDiS/src/Lagrange.cc index a1a8d835ef11c627522f207624086d361e202e9e..9fbacbd55e857ac4c1bd535e5cf946e6393eb9ff 100644 --- a/AMDiS/src/Lagrange.cc +++ b/AMDiS/src/Lagrange.cc @@ -792,7 +792,7 @@ namespace AMDiS { DimVec<int> parts(dim, NO_INIT); for (int i = 0; i < dim + 1; i++) - parts[i] = Global::getGeo(INDEX_OF_DIM(i, dim),dim); + parts[i] = Global::getGeo(INDEX_OF_DIM(i, dim), dim); int index = 0; diff --git a/AMDiS/src/Mesh.cc b/AMDiS/src/Mesh.cc index b781ab73568ea41a8a9d8d5d9c36ba0fd096e2be..0e31ba4152ff575658143ac60956096edd4065ae 100644 --- a/AMDiS/src/Mesh.cc +++ b/AMDiS/src/Mesh.cc @@ -617,25 +617,23 @@ namespace AMDiS { static const MacroElement *mel = NULL; DimVec<double> lambda(dim, NO_INIT); ElInfo *mel_info = NULL; - int i, k; - bool inside; mel_info = createNewElInfo(); if (start_mel != NULL) mel = start_mel; else - if((mel == NULL)||(mel->getElement()->getMesh() != this)) + if ((mel == NULL) || (mel->getElement()->getMesh() != this)) mel = *(macroElements.begin()); mel_info->setFillFlag(Mesh::FILL_COORDS); - g_xy = &xy; + g_xy = &xy; g_xy0 = xy0; - g_sp = sp; + g_sp = sp; mel_info->fillMacroInfo(mel); - + int k; while ((k = mel_info->worldToCoord(xy, &lambda)) >= 0) { if (mel->getNeighbour(k)) { mel = mel->getNeighbour(k); @@ -646,8 +644,10 @@ namespace AMDiS { } /* now, descend in tree to find leaf element at point */ - inside = findElementAtPointRecursive(mel_info, lambda, k, el_info); - for (i=0; i<=dim; i++) bary[i] = final_lambda[i]; + bool inside = findElementAtPointRecursive(mel_info, lambda, k, el_info); + for (int i = 0; i <= dim; i++) { + bary[i] = final_lambda[i]; + } DELETE mel_info; @@ -655,18 +655,14 @@ namespace AMDiS { } bool Mesh::findElementAtPoint(const WorldVector<double>& xy, - Element **elp, - DimVec<double>& bary, + Element **elp, + DimVec<double>& bary, const MacroElement *start_mel, - const WorldVector<double> *xy0, - double *sp) + const WorldVector<double> *xy0, + double *sp) { - ElInfo *el_info = NULL; - int val; - - el_info = createNewElInfo(); - - val = findElInfoAtPoint(xy, el_info, bary, start_mel, xy0, sp); + ElInfo *el_info = createNewElInfo(); + int val = findElInfoAtPoint(xy, el_info, bary, start_mel, xy0, sp); *elp = el_info->getElement(); @@ -677,48 +673,47 @@ namespace AMDiS { - bool Mesh::findElementAtPointRecursive(ElInfo *el_info, + bool Mesh::findElementAtPointRecursive(ElInfo *el_info, const DimVec<double>& lambda, - int outside, + int outside, ElInfo* final_el_info) { - FUNCNAME("Mesh::findElementAtPointRecursive"); + FUNCNAME("Mesh::findElementAtPointRecursive()"); Element *el = el_info->getElement(); - ElInfo *c_el_info = NULL; DimVec<double> c_lambda(dim, NO_INIT); - int i, inside; - int ichild, c_outside; + int inside; + int ichild, c_outside; if (el->isLeaf()) { *final_el_info = *el_info; if (outside < 0) { - for (i=0; i<=dim; i++) final_lambda[i] = lambda[i]; + for (int i = 0; i <= dim; i++) { + final_lambda[i] = lambda[i]; + } + return(true); - } - else - { /* outside */ - if (g_xy0) - { /* find boundary point of [xy0, xy] */ - double s; - el_info->worldToCoord(*(g_xy0), &c_lambda); - s = lambda[outside] / (lambda[outside] - c_lambda[outside]); - for (i=0; i<=dim; i++) - { - final_lambda[i] = s * c_lambda[i] + (1.0-s) * lambda[i]; - } - if (g_sp) *(g_sp) = s; - if(dim == 3) - MSG("outside finest level on el %d: s=%.3e\n", el->getIndex(), s); + } else { /* outside */ + if (g_xy0) { /* find boundary point of [xy0, xy] */ + el_info->worldToCoord(*(g_xy0), &c_lambda); + double s = lambda[outside] / (lambda[outside] - c_lambda[outside]); + for (int i = 0; i <= dim; i++) { + final_lambda[i] = s * c_lambda[i] + (1.0-s) * lambda[i]; + } + if (g_sp) { + *(g_sp) = s; + } + if (dim == 3) + MSG("outside finest level on el %d: s=%.3e\n", el->getIndex(), s); - return(false); /* ??? */ - } - else return(false); + return(false); /* ??? */ } + else return(false); + } } - c_el_info = createNewElInfo(); + ElInfo *c_el_info = createNewElInfo(); - if(dim == 1) { + if (dim == 1) { if (lambda[0] >= lambda[1]) { c_el_info->fillElInfo(0, el_info); if (outside >= 0) { @@ -740,7 +735,7 @@ namespace AMDiS { } } /* DIM == 1 */ - if(dim == 2) { + if (dim == 2) { if (lambda[0] >= lambda[1]) { c_el_info->fillElInfo(0, el_info); if (el->isNewCoordSet()) { @@ -768,7 +763,7 @@ namespace AMDiS { } } /* DIM == 2 */ - if(dim == 3) { + if (dim == 3) { if (el->isNewCoordSet()) { if (lambda[0] >= lambda[1]) ichild = 0; @@ -792,7 +787,9 @@ namespace AMDiS { c_lambda2[3]); if ((c_outside2 < 0) || (c_lambda2[c_outside2] > c_lambda[c_outside])) { - for (i=0; i<=dim; i++) c_lambda[i] = c_lambda2[i]; + for (int i = 0; i <= dim; i++) { + c_lambda[i] = c_lambda2[i]; + } c_outside = c_outside2; *c_el_info = *c_el_info2; ichild = 1 - ichild; diff --git a/AMDiS/src/OpenMP.h b/AMDiS/src/OpenMP.h index 455e1838dd43c657830ea70288c305f757586969..63cd147377e63fbf18609f35e4adca9e057c4011 100644 --- a/AMDiS/src/OpenMP.h +++ b/AMDiS/src/OpenMP.h @@ -5,8 +5,12 @@ #include <omp.h> +const bool amdisHaveOpenMP = true; + #else +const bool amdisHaveOpenMP = false; + inline int omp_get_max_threads() { return 1; } diff --git a/AMDiS/src/Preconditioner.h b/AMDiS/src/Preconditioner.h index ccea6bd7aeed0ee5205d6273a6086d8171f31c1d..21dc117e7d3967a09637c5d844fc59e1618c1a60 100644 --- a/AMDiS/src/Preconditioner.h +++ b/AMDiS/src/Preconditioner.h @@ -240,7 +240,9 @@ namespace AMDiS { int i; int size = scalPrecons.getSize(); +#ifdef _OPENMP #pragma omp parallel for +#endif for (i = 0; i < size; i++) { scalPrecons[i]->init(); } @@ -252,8 +254,9 @@ namespace AMDiS { virtual void precon(SystemVector *x) { int i; int size = scalPrecons.getSize(); - +#ifdef _OPENMP #pragma omp parallel for +#endif for (i = 0; i < size; i++) { scalPrecons[i]->precon(x->getDOFVector(i)); } @@ -265,8 +268,9 @@ namespace AMDiS { virtual void exit() { int i; int size = scalPrecons.getSize(); - +#ifdef _OPENMP #pragma omp parallel for +#endif for (i = 0; i < size; i++) { scalPrecons[i]->exit(); } diff --git a/AMDiS/src/Tetrahedron.h b/AMDiS/src/Tetrahedron.h index c638717f668b5f39eae9ebcc158621263533b7dc..ceb031f37a4bc01fad85d4750f4b551b7cba4b81 100644 --- a/AMDiS/src/Tetrahedron.h +++ b/AMDiS/src/Tetrahedron.h @@ -48,7 +48,9 @@ namespace AMDiS { /** \brief * implements Element::clone */ - inline Element *clone() { return NEW Tetrahedron(mesh); }; + inline Element *clone() { + return NEW Tetrahedron(mesh); + }; /** \brief * implements Element::getVertexOfEdge @@ -103,13 +105,6 @@ namespace AMDiS { } }; - /** \brief - * implements Element::sideOfChild - */ - //int sideOfChild(int child, - // int side, bool *isVirtual, - // unsigned char elementType=0); - /** \brief * implements Element::hasSide */ @@ -125,23 +120,31 @@ namespace AMDiS { * implements Element::isLine. Returns false because this element is a * Tetrahedron */ - inline bool isLine() const { return false; }; + inline bool isLine() const { + return false; + }; /** \brief * implements Element::isTriangle. Returns false because this element is a * Tetrahedron */ - inline bool isTriangle() const { return false; }; + inline bool isTriangle() const { + return false; + }; /** \brief * implements Element::isTetrahedron. Returns true because this element is a * Tetrahedron */ - inline bool isTetrahedron() const { return true; }; + inline bool isTetrahedron() const { + return true; + }; // ===== Serializable implementation ===== - ::std::string getTypeName() const { return "Tetrahedron"; }; + ::std::string getTypeName() const { + return "Tetrahedron"; + }; public: /** \brief