Skip to content
Snippets Groups Projects
Commit ea26c5eb authored by Praetorius, Simon's avatar Praetorius, Simon
Browse files

correct parallel evalAtPoint functions if point is found on multiple partitions

parent 8aa43b17
No related branches found
No related tags found
2 merge requests!28Release/v1.2,!27correct parallel evalAtPoint functions if point is found on multiple partitions
......@@ -124,10 +124,10 @@ namespace AMDiS {
{
this->name = n;
this->feSpace = f;
if (this->feSpace && this->feSpace->getAdmin())
(this->feSpace->getAdmin())->addDOFIndexed(this);
this->boundaryManager = new BoundaryManager(f);
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
if (addToSynch && Parallel::MeshDistributor::globalMeshDistributor != NULL)
......@@ -142,7 +142,7 @@ namespace AMDiS {
if ( Parallel::MeshDistributor::globalMeshDistributor != NULL)
Parallel::MeshDistributor::globalMeshDistributor->removeInterchangeVector(this);
#endif
#ifdef _OPENMP
#pragma omp critical
#endif
......@@ -1090,7 +1090,7 @@ namespace AMDiS {
ElInfo *elInfo = mesh->createNewElInfo();
double value = 0.0;
bool inside = false;
unsigned short inside = 0;
if (oldElInfo && oldElInfo->getMacroElement()) {
inside = mesh->findElInfoAtPoint(p, elInfo, lambda, oldElInfo->getMacroElement(), NULL, NULL);
......@@ -1101,28 +1101,24 @@ namespace AMDiS {
if (oldElInfo)
oldElInfo = elInfo;
if (inside) {
if (inside > 0) {
basFcts->getLocalIndices(elInfo->getElement(), this->feSpace->getAdmin(), localIndices);
ElementVector uh(nBasFcts);
for (int i = 0; i < nBasFcts; i++)
uh[i] = operator[](localIndices[i]);
value = basFcts->evalUh(lambda, uh);
} else {
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
value = 0.0;
#else
ERROR_EXIT("Can not eval DOFVector at point p, because point is outside geometry.");
#endif
}
if (oldElInfo == NULL)
delete elInfo;
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
Parallel::mpi::globalAdd(value);
Parallel::mpi::globalAdd(inside);
#endif
return value;
TEST_EXIT(inside > 0)("Can not eval DOFVector at point p, because point is outside geometry.");
return value/inside;
}
......@@ -1140,7 +1136,7 @@ namespace AMDiS {
DimVec<double> lambda(dim, NO_INIT);
ElInfo *elInfo = mesh->createNewElInfo();
WorldVector<double> value(DEFAULT_VALUE, 0.0);
bool inside = false;
unsigned short inside = 0;
if (oldElInfo && oldElInfo->getMacroElement()) {
inside = mesh->findElInfoAtPoint(p, elInfo, lambda, oldElInfo->getMacroElement(), NULL, NULL);
......@@ -1151,19 +1147,24 @@ namespace AMDiS {
if (oldElInfo)
oldElInfo = elInfo;
if (inside) {
if (inside > 0) {
basFcts->getLocalIndices(elInfo->getElement(), this->feSpace->getAdmin(), localIndices);
mtl::dense_vector<WorldVector<double> > uh(nBasFcts);
for (int i = 0; i < nBasFcts; i++)
uh[i] = operator[](localIndices[i]);
value = basFcts->evalUh(lambda, uh);
} else {
ERROR_EXIT("Can not eval DOFVector at point p, because point is outside geometry.");
}
if (oldElInfo == NULL)
delete elInfo;
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
Parallel::mpi::globalAdd(value);
Parallel::mpi::globalAdd(inside);
#endif
TEST_EXIT(inside > 0)("Can not eval DOFVector at point p, because point is outside geometry.");
value *= 1.0/inside;
return value;
}
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment