Commit 1500f6ae authored by Thomas Witkowski's avatar Thomas Witkowski

* Dies und jenes

* OpenMP-Versuche aus Vorkonditionierern und Loesern geloescht
parent afc60495
......@@ -331,7 +331,7 @@ link_all_deplibs=unknown
sys_lib_search_path_spec=`echo " /usr/lib/gcc/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../../i386-redhat-linux/lib/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../../i386-redhat-linux/lib/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../ /lib/i386-redhat-linux/4.1.2/ /lib/ /usr/lib/i386-redhat-linux/4.1.2/ /usr/lib/" | $SED -e "s@${gcc_dir}@\${gcc_dir}@g;s@${gcc_ver}@\${gcc_ver}@g"`
# Run-time system search path for libraries
sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/qt-3.3/lib "
sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/octave-2.9.9 /usr/lib/qt-3.3/lib "
# Fix the shell variable $srcfile for the compiler.
fix_srcfile_path=""
......@@ -7550,7 +7550,7 @@ link_all_deplibs=unknown
sys_lib_search_path_spec=`echo " /usr/lib/gcc/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../../i386-redhat-linux/lib/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../../i386-redhat-linux/lib/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../ /lib/i386-redhat-linux/4.1.2/ /lib/ /usr/lib/i386-redhat-linux/4.1.2/ /usr/lib/" | $SED -e "s@${gcc_dir}@\${gcc_dir}@g;s@${gcc_ver}@\${gcc_ver}@g"`
# Run-time system search path for libraries
sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/qt-3.3/lib "
sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/octave-2.9.9 /usr/lib/qt-3.3/lib "
# Fix the shell variable $srcfile for the compiler.
fix_srcfile_path=""
......@@ -7861,7 +7861,7 @@ link_all_deplibs=unknown
sys_lib_search_path_spec=`echo " /usr/lib/gcc/i386-redhat-linux/3.4.6/ /usr/lib/gcc/i386-redhat-linux/3.4.6/ /usr/lib/gcc/i386-redhat-linux/3.4.6/../../../../i386-redhat-linux/lib/i386-redhat-linux/3.4.6/ /usr/lib/gcc/i386-redhat-linux/3.4.6/../../../../i386-redhat-linux/lib/ /usr/lib/gcc/i386-redhat-linux/3.4.6/../../../i386-redhat-linux/3.4.6/ /usr/lib/gcc/i386-redhat-linux/3.4.6/../../../ /lib/i386-redhat-linux/3.4.6/ /lib/ /usr/lib/i386-redhat-linux/3.4.6/ /usr/lib/" | $SED -e "s@${gcc_dir}@\${gcc_dir}@g;s@${gcc_ver}@\${gcc_ver}@g"`
# Run-time system search path for libraries
sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/qt-3.3/lib "
sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/octave-2.9.9 /usr/lib/qt-3.3/lib "
# Fix the shell variable $srcfile for the compiler.
fix_srcfile_path=""
......
......@@ -497,6 +497,7 @@ namespace AMDiS {
int numPoints = quadrature->getNumPoints();
double *c = GET_MEMORY(double, numPoints);
for (int iq = 0; iq < numPoints; iq++) {
c[iq] = 0.0;
}
......@@ -525,7 +526,7 @@ namespace AMDiS {
}
}
}
} else { /* non symmetric assembling */
} else { // non symmetric assembling
for (int iq = 0; iq < numPoints; iq++) {
c[iq] *= elInfo->getDet();
......@@ -542,6 +543,7 @@ namespace AMDiS {
}
}
}
FREE_MEMORY(c, double, numPoints);
FREE_MEMORY(phival, double, nCol);
}
......@@ -1394,10 +1396,11 @@ namespace AMDiS {
Element *el = elInfo->getElement();
checkForNewTraverse();
checkQuadratures();
if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) {
initElement(elInfo);
}
......@@ -1410,8 +1413,6 @@ namespace AMDiS {
} else {
if (rememberElMat) {
axpy(factor, *elementMatrix, *userMat);
//*userMat += *elementMatrix * factor;
//operat->addElementMatrix(elementMatrix, userMat, factor);
return;
}
}
......@@ -1427,11 +1428,9 @@ namespace AMDiS {
if (zeroOrderAssembler)
zeroOrderAssembler->calculateElementMatrix(elInfo, mat);
if(rememberElMat && userMat) {
if (rememberElMat && userMat) {
axpy(factor, *elementMatrix, *userMat);
//*userMat += *elementMatrix * factor;
//operat->addElementMatrix(elementMatrix, userMat, factor);
}
}
}
void Assembler::calculateElementVector(const ElInfo *elInfo,
......@@ -1597,27 +1596,40 @@ namespace AMDiS {
ElementMatrix *Assembler::initElementMatrix(ElementMatrix *elMat,
const ElInfo *elInfo)
{
if(!elMat) {
if (!elMat) {
elMat = NEW ElementMatrix(nRow, nCol);
}
elMat->set(0.0);
DOFAdmin *rowAdmin = rowFESpace->getAdmin();
DOFAdmin *colAdmin = colFESpace->getAdmin();
Element *element = elInfo->getElement();
/*
elMat->rowIndices =
rowFESpace->getBasisFcts()->getLocalIndices(element,
rowAdmin,
NULL);
*/
/*
elMat->colIndices =
colFESpace->getBasisFcts()->getLocalIndices(element,
colAdmin,
NULL);
*/
rowFESpace->getBasisFcts()->getLocalIndicesVec(element,
rowAdmin,
&(elMat->rowIndices));
colFESpace->getBasisFcts()->getLocalIndicesVec(element,
colAdmin,
&(elMat->colIndices));
return elMat;
}
......
......@@ -29,6 +29,7 @@
#include <string>
#include "Global.h"
#include "Boundary.h"
#include "MatrixVector.h"
namespace AMDiS {
......@@ -318,6 +319,14 @@ namespace AMDiS {
return NULL;
};
/** \brief
* Returns local dof indices of the element for the given fe space.
*/
virtual void getLocalIndicesVec(const Element*,
const DOFAdmin*,
Vector<DegreeOfFreedom>*) const
{};
/** \brief
* Evaluates elements value at barycentric coordinates lambda with local
......
......@@ -450,15 +450,14 @@ namespace AMDiS {
FUNCNAME("DOFMatrix::assemble()");
if (!op && operators.size() == 0) {
//WARNING("no operator\n");
return NULL;
}
Operator *operat = op ? op : operators[0];
elementMatrix =
operat->getAssembler()->initElementMatrix(elementMatrix, elInfo);
if (op) {
op->getElementMatrix(elInfo, elementMatrix);
} else {
......@@ -472,10 +471,11 @@ namespace AMDiS {
elementMatrix,
*factorIt ? **factorIt : 1.0);
}
}
addElementMatrix(factor, *elementMatrix, bound);
return elementMatrix;
}
......
......@@ -1207,61 +1207,29 @@ namespace AMDiS {
("a.size = %d too small: admin->sizeUsed = %d\n",
a.getSize(), rowAdmin->getUsedSize());
// This is the old implementation of the mv-multiplication. It have been changed
// because of the OpenMP-parallelization:
// typename DOFVector<T>::Iterator vecIterator(dynamic_cast<DOFIndexed<T>*>(&result), USED_DOFS);
// DOFMatrix::Iterator rowIterator(const_cast<DOFMatrix*>(&a), USED_DOFS);
// for(vecIterator.reset(), rowIterator.reset();
// !rowIterator.end();
// ++rowIterator, ++vecIterator) {
// sum = 0;
// if(!add) *vecIterator = 0.0;
// for(::std::vector<MatEntry>::iterator colIterator = rowIterator->begin();
// colIterator != rowIterator->end();
// colIterator++) {
// jcol = colIterator->col;
// if (jcol >= 0) { // entry used?
// sum += (static_cast<double>(colIterator->entry)) * x[jcol];
// } else {
// if (jcol == DOFMatrix::NO_MORE_ENTRIES)
// break;
// }
// }
// *vecIterator += sum;
// }
int i;
int maxI = result.getSize();
#ifdef _OPENMP
#pragma omp parallel for schedule(dynamic, 25000) default(shared) private(i)
#endif
for (i = 0; i < maxI; i++) {
if (rowAdmin->isDOFFree(i)) {
continue;
}
T sum = 0;
const ::std::vector<MatEntry> *v = &a[i];
for (int j = 0; j < static_cast<int>(v->size()); j++) {
const MatEntry *m = &((*v)[j]);
int jcol = m->col;
if (jcol >= 0) { // entry used?
sum += (static_cast<double>(m->entry)) * x[jcol];
} else {
if (jcol == DOFMatrix::NO_MORE_ENTRIES)
break;
}
}
if (add) {
result[i] += sum;
} else {
result[i] = sum;
}
typename DOFVector<T>::Iterator vecIterator(dynamic_cast<DOFIndexed<T>*>(&result), USED_DOFS);
DOFMatrix::Iterator rowIterator(const_cast<DOFMatrix*>(&a), USED_DOFS);
for(vecIterator.reset(), rowIterator.reset();
!rowIterator.end();
++rowIterator, ++vecIterator) {
double sum = 0;
if (!add)
*vecIterator = 0.0;
for(::std::vector<MatEntry>::iterator colIterator = rowIterator->begin();
colIterator != rowIterator->end();
colIterator++) {
int jcol = colIterator->col;
if (jcol >= 0) { // entry used?
sum += (static_cast<double>(colIterator->entry)) * x[jcol];
} else {
if (jcol == DOFMatrix::NO_MORE_ENTRIES)
break;
}
}
*vecIterator += sum;
}
} else if (transpose == Transpose) {
......
......@@ -235,10 +235,6 @@ namespace AMDiS {
el = elInfo->getElement();
// if (el->getIndex() == 66065) {
// ::std::cout << "UNSER ELEMENT!" << ::std::endl;
// }
double det = elInfo->getDet();
const DimVec<WorldVector<double> > &Lambda = elInfo->getGrdLambda();
......@@ -246,11 +242,11 @@ namespace AMDiS {
double h2 = h2_from_det(det, dim);
for(iq = 0; iq < numPoints; iq++) {
for (iq = 0; iq < numPoints; iq++) {
riq[iq] = 0.0;
}
for(system = 0; system < numSystems; system++) {
for (system = 0; system < numSystems; system++) {
if(matrix[system] == NULL) continue;
......
......@@ -1004,7 +1004,7 @@ namespace AMDiS {
GeoIndex posIndex;
DegreeOfFreedom* result;
if (indices) {
result = indices;
} else {
......@@ -1014,7 +1014,7 @@ namespace AMDiS {
localVecSize = nBasFcts;
result = localVec;
}
for (pos = 0, j = 0; pos <= dim; pos++) {
posIndex = INDEX_OF_DIM(pos, dim);
n0 = admin->getNumberOfPreDOFs(posIndex);
......@@ -1031,10 +1031,21 @@ namespace AMDiS {
}
}
}
return result;
}
void Lagrange::getLocalIndicesVec(const Element* el,
const DOFAdmin *admin,
Vector<DegreeOfFreedom> *indices) const
{
if (indices->getSize() < nBasFcts) {
indices->resize(nBasFcts);
}
getLocalIndices(el, admin, &((*indices)[0]));
}
void Lagrange::l2ScpFctBas(Quadrature *q,
AbstractFunction<WorldVector<double>, WorldVector<double> >* f,
DOFVector<WorldVector<double> >* fh)
......
......@@ -152,6 +152,11 @@ namespace AMDiS {
const DOFAdmin*,
DegreeOfFreedom*) const;
void getLocalIndicesVec(const Element*,
const DOFAdmin*,
Vector<DegreeOfFreedom>*) const;
/** \brief
* Implements BasisFunction::l2ScpFctBas
*/
......
......@@ -237,11 +237,7 @@ namespace AMDiS {
* Initialisation of the preconditioner
*/
virtual void init() {
int i;
#ifdef _OPENMP
#pragma omp parallel for schedule(static, 1) default(shared) private(i)
#endif
for (i = 0; i < scalPrecons.getSize(); i++) {
for (int i = 0; i < scalPrecons.getSize(); i++) {
scalPrecons[i]->init();
}
};
......@@ -250,11 +246,7 @@ namespace AMDiS {
* Preconditioning method
*/
virtual void precon(SystemVector *x) {
int i;
#ifdef _OPENMP
#pragma omp parallel for schedule(static, 1) default(shared) private(i)
#endif
for (i = 0; i < scalPrecons.getSize(); i++) {
for (int i = 0; i < scalPrecons.getSize(); i++) {
scalPrecons[i]->precon(x->getDOFVector(i));
}
};
......@@ -263,11 +255,7 @@ namespace AMDiS {
* Frees needed memory.
*/
virtual void exit() {
int i;
#ifdef _OPENMP
#pragma omp parallel for schedule(static, 1) default(shared) private(i)
#endif
for (i = 0; i < scalPrecons.getSize(); i++) {
for (int i = 0; i < scalPrecons.getSize(); i++) {
scalPrecons[i]->exit();
}
};
......
......@@ -679,10 +679,6 @@ namespace AMDiS {
}
}
TraverseStack stack;
DualTraverse dualStack;
ElInfo *elInfo;
for (int i = 0; i < numComponents_; i++) {
for (int j = 0; j < numComponents_; j++) {
// Only if this variable is true, the current matrix will be assembled.
......@@ -708,19 +704,22 @@ namespace AMDiS {
continue;
}
if (componentMeshes_[i] != componentMeshes_[j]) {
ERROR_EXIT("not yet\n");
} else {
if (assembleMatrix && matrix->getBoundaryManager())
matrix->getBoundaryManager()->initMatrix(matrix);
elInfo = stack.traverseFirst(componentMeshes_[i], -1, assembleFlag);
TraverseStack stack;
ElInfo *elInfo = stack.traverseFirst(componentMeshes_[i], -1, assembleFlag);
while (elInfo) {
const BoundaryType *bound =
useGetBound_ ?
componentSpaces_[i]->getBasisFcts()->getBound(elInfo, NULL) :
NULL;
if (assembleMatrix) {
matrix->assemble(1.0, elInfo, bound);
if (matrix->getBoundaryManager()) {
......@@ -732,12 +731,14 @@ namespace AMDiS {
if (i == j) {
rhs_->getDOFVector(i)->assemble(1.0, elInfo, bound);
}
elInfo = stack.traverseNext(elInfo);
}
if (assembleMatrix && matrix->getBoundaryManager())
matrix->getBoundaryManager()->exitMatrix(matrix);
matrix->getBoundaryManager()->exitMatrix(matrix);
}
assembledMatrix_[i][j] = true;
}
......@@ -748,7 +749,8 @@ namespace AMDiS {
if (solution_->getDOFVector(i)->getBoundaryManager())
solution_->getDOFVector(i)->getBoundaryManager()->initVector(solution_->getDOFVector(i));
elInfo = stack.traverseFirst(componentMeshes_[i], -1, assembleFlag);
TraverseStack stack;
ElInfo *elInfo = stack.traverseFirst(componentMeshes_[i], -1, assembleFlag);
while (elInfo) {
if(rhs_->getDOFVector(i)->getBoundaryManager())
rhs_->getDOFVector(i)->getBoundaryManager()->
......@@ -761,28 +763,11 @@ namespace AMDiS {
if (rhs_->getDOFVector(i)->getBoundaryManager())
rhs_->getDOFVector(i)->getBoundaryManager()->exitVector(rhs_->getDOFVector(i));
if (solution_->getDOFVector(i)->getBoundaryManager())
solution_->getDOFVector(i)->getBoundaryManager()->exitVector(solution_->getDOFVector(i));
solution_->getDOFVector(i)->getBoundaryManager()->exitVector(solution_->getDOFVector(i));
}
for (int i = 0; i < numComponents_; i++) {
// ::std::cout << "Vector: " << i << ::std::endl;
// rhs_->getDOFVector(i)->print();
}
INFO(info_, 8)("buildAfterCoarsen needed %.5f seconds\n",
TIME_USED(first,clock()));
// int memsize = 0;
// for (int i = 0; i < numComponents_; i++) {
// for (int j = 0; j < numComponents_; j++) {
// DOFMatrix *matrix = (*systemMatrix_)[i][j];
// memsize += (*matrix).memsize();
// }
// }
// ::std::cout << "SIZE: " << memsize << ::std::endl;
TIME_USED(first, clock()));
}
void ProblemVec::writeFiles(AdaptInfo *adaptInfo, bool force)
......
......@@ -474,14 +474,11 @@ namespace AMDiS {
TEST_EXIT(size == matrix.getNumRows())("incompatible sizes\n");
TEST_EXIT(size == matrix.getNumCols())("incompatible sizes\n");
int i, j;
#ifdef _OPENMP
#pragma omp parallel for schedule(static, 1) default(shared) private(i,j)
#endif
for (i = 0; i < size; i++) {
if (!add) result.getDOFVector(i)->set(0.0);
for (j = 0; j < size; j++) {
for (int i = 0; i < size; i++) {
if (!add)
result.getDOFVector(i)->set(0.0);
for (int j = 0; j < size; j++) {
if (matrix[i][j]) {
mv<double>(NoTranspose,
*(matrix[i][j]),
......@@ -500,12 +497,10 @@ namespace AMDiS {
{
TEST_EXIT(x.getNumVectors() == y.getNumVectors())
("invalid size\n");
int i, size = x.getNumVectors();
#ifdef _OPENMP
#pragma omp parallel for schedule(static, 1) default(shared) private(i)
#endif
for (i = 0; i < size; i++) {
int size = x.getNumVectors();
for (int i = 0; i < size; i++) {
axpy(a, *(x.getDOFVector(i)), *(y.getDOFVector(i)));
}
};
......@@ -517,12 +512,9 @@ namespace AMDiS {
{
TEST_EXIT(x.getNumVectors() == y.getNumVectors())
("invalid size\n");
int i, size = x.getNumVectors();
int size = x.getNumVectors();
#ifdef _OPENMP
#pragma omp parallel for schedule(static, 1) default(shared) private(i)
#endif
for (i = 0; i < size; i++) {
for (int i = 0; i < size; i++) {
xpay(a, *(x.getDOFVector(i)), *(y.getDOFVector(i)));
}
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment