Homec4science
Diffusion LibMultiScale (master) Browse

libmultiscale/src/coupling/point_association.cc
41eb1f366000master

Last Change
Raw File

point_association.cc
View Options

/**
* @file point_association.cc
*
* @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
* @author Till Junge <till.junge@epfl.ch>
*
* @date Thu Jul 24 14:21:58 2014
*
* @brief Coupling tool to associate point-like DOFs
*
* @section LICENSE
*
* Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* LibMultiScale is free software: you can redistribute it and/or modify it
* under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* LibMultiScale is distributed in the hope that it will be useful, but
* WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with LibMultiScale. If not, see <http://www.gnu.org/licenses/>.
*
*/
//#define TIMER
/* -------------------------------------------------------------------------- */
#include "point_association.hh"
#include "filter_geometry.hh"
#include "geometry_manager.hh"
#include "lib_bridging.hh"
#include "lib_continuum.hh"
#include "lib_md.hh"
#include "lm_common.hh"
#include "reference_manager_interface.hh"
/* -------------------------------------------------------------------------- */
__BEGIN_LIBMULTISCALE__
/* -------------------------------------------------------------------------- */
PointAssociation::PointAssociation(const std::string &name)
: LMObject(name), DofAssociation(name), pointsAlist("pointsAlist:" + name),
pointsBlist("pointsBlist:" + name), unmatchedA("unmatchedA:" + name),
unmatchedB("unmatchedB:" + name), assoc_found(0), tolerance(1e-5),
path(""), weight_path(""), weights("weights:" + this->getID()),
weightsA("weightsA:" + this->getID()),
weightsB("weightsB:" + this->getID()),
true_disp("ComputeTrueDisplacement:" + this->getID()),
grid(createOutput("grid")) {}
/* -------------------------------------------------------------------------- */
PointAssociation::~PointAssociation() {}
/* -------------------------------------------------------------------------- */
void PointAssociation::fillGrid() {
LM_TOIMPLEMENT;
// constexpr UInt Dim = ContainerA::Dim;
// create my grid
// this->grid.reset();
// Geometry &geom = *GeometryManager::getManager().getGeometry(this->geomID);
// Cube cube = geom.getBoundingBox();
// LM_TOIMPLEMENT;
// // this->grid = new SpatialGridPoint<UInt, Dim>(cube, this->grid_division);
// UInt index = 0;
// for (auto &&point : this->unmatchedB) {
// auto position = point.position0();
// // I'ts important to use addAtom here instead of addElement, because here
// // we're sensitive to atoms on the edge of grid_boxes
// this->grid->addAtom(index, position);
// ++index;
// }
}
/* -------------------------------------------------------------------------- */
template <UInt Dim, typename V1, typename V2>
inline bool coincidence(V1 &pointA, V2 &pointB, const Real &tolerance) {
Real distance = 0;
for (UInt dir = 0; dir < Dim; ++dir) {
distance += fabs(pointA[dir] - pointB[dir]);
}
return distance < tolerance;
}
/* -------------------------------------------------------------------------- */
template <UInt Dim, typename V1, typename V2>
inline Real temp_coincidence(V1 &pointA, V2 &pointB) {
Real distance = 0;
for (UInt dir = 0; dir < Dim; ++dir) {
distance += fabs(pointA[dir] - pointB[dir]);
}
return distance;
}
/* -------------------------------------------------------------------------- */
void PointAssociation::associate() {
LM_TOIMPLEMENT;
// this->pointToElement.assign(this->local_points, lm_uint_max);
// constexpr UInt Dim = ContainerA::Dim;
// LM_ASSERT(this->positions.size() == this->local_points * Dim,
// "inconsistency " << this->positions.size()
// << " != " << this->local_points * Dim);
// // Loop over all positions of A world
// bool found;
// Real min_dist;
// UInt not_found_counter = 0;
// Real *posA_ptr = &this->positions[0];
// for (UInt indexA = 0; indexA < this->local_points;
// ++indexA, posA_ptr += Dim) {
// // extract the points of pointsBlist which are candidates for association
// VectorView<Dim> positionA(posA_ptr);
// LM_TOIMPLEMENT;
// // auto &&subsetPointsB = this->grid->findSet(positionA);
// found = false;
// min_dist = lm_real_max;
// // loop the candidates
// LM_TOIMPLEMENT;
// // for (auto &&indexB : subsetPointsB) {
// // auto &&positionB = unmatchedB.get(indexB).position0();
// // Real dist = temp_coincidence<Dim>(positionA, positionB);
// // if (dist < min_dist)
// // min_dist = dist;
// // if (coincidence<Dim>(positionA, positionB, this->tolerance)) {
// // this->assoc[indexA] = indexB;
// // ++this->assoc_found;
// // found = true;
// // if (indexA < 10)
// // DUMP("associate point " << positionA << ", " << indexB
// // << " with position " << indexA,
// // DBG_DETAIL);
// // break;
// // }
// // }
// if (!found) {
// ++not_found_counter;
// {
// std::stringstream sstr;
// sstr << not_found_counter << ": couldn't associate point (";
// for (UInt i = 0; i < Dim - 1; ++i)
// sstr << positionA[i] << ", ";
// sstr << positionA[Dim - 1] << "), "
// << "the minimal distance was " << min_dist;
// DUMP(sstr.str(), DBG_INFO);
// }
// }
// }
}
/* -------------------------------------------------------------------------- */
void PointAssociation::clearAll() {
this->assoc_found = 0;
this->point2Point.clear();
LM_TOIMPLEMENT;
// this->grid.reset();
// pointsAlist->clear();
// pointsBlist->clear();
}
/* -------------------------------------------------------------------------- */
void PointAssociation::filterPointListAForUnmatched() {
LM_TOIMPLEMENT;
// this->pointsAlist->clear();
// for (UInt i = 0; i < this->pointToElement.size(); ++i) {
// if (this->pointToElement[i] != UINT_MAX) {
// DUMP("for container, atom " << i << " becomes "
// << this->pointsAlist->size(),
// DBG_ALL);
// this->pointsAlist->push_back(this->unmatchedA->get(i));
// } else {
// DUMP("atom filtered " << i << " " << this->unmatchedA->get(i),
// DBG_DETAIL);
// }
// }
// this->pointsAlist->copyContainerInfo(unmatchedA);
LM_TOIMPLEMENT;
// this->pointsAlist->copyReleaseInfo(unmatchedA);
}
/* -------------------------------------------------------------------------- */
void PointAssociation::filterPointListBForUnmatched() {
LM_TOIMPLEMENT;
// this->pointsBlist->clear();
// UInt counter = 0;
// for (UInt i = 0; i < this->pointToElement.size(); ++i) {
// if (this->pointToElement[i] != UINT_MAX) {
// DUMP("for container, atom " << i << " becomes "
// << this->pointsBlist->size(),
// DBG_ALL);
// this->pointsBlist->push_back(unmatchedB->get(this->pointToElement[i]));
// this->pointToElement[i] = counter;
// ++counter;
// } else {
// DUMP("atom filtered " << i, DBG_DETAIL);
// }
// }
// this->pointsBlist->getBoundingBox() = unmatchedB->getBoundingBox();
}
/* -------------------------------------------------------------------------- */
void PointAssociation::checkPositionCoherency() {
LM_TOIMPLEMENT;
// try {
// std::stringstream sstr;
// sstr << this->getID() << "-duovector-" << current_step << "-" <<
// current_stage;
// this->getDuoVector().print(sstr.str());
// }
// catch (LibMultiScaleException & e){
// }
// extract the position of local atoms
// if (this->in_group_A)
// this->buildPositions(this->pointsAlist.getOutput());
// constexpr UInt Dim = ContainerA::Dim;
// // transport the positions to B side
// this->distributeVectorA2B("temp_positions", this->positions);
// if (this->in_group_B) {
// UInt point_idx = 0;
// for (auto &&point : this->pointsBlist) {
// auto &&local_pos = point.position0();
// std::stringstream sstr;
// for (UInt i = 0; i < Dim; ++i)
// sstr << this->positions[Dim * point_idx + i] << " ";
// sstr << " and ";
// for (UInt i = 0; i < Dim; ++i)
// sstr << local_pos[i] << " ";
// std::string _tmp = sstr.str();
// DUMP("for point " << point_idx << " compare " << _tmp, DBG_DETAIL);
// for (UInt k = 0; k < Dim; ++k) {
// if (fabs(this->positions[Dim * point_idx + k] - local_pos[k]) > 1e-5)
// {
// LM_FATAL(this->getID()
// << ": inconsistency in the communication scheme "
// << this->positions[Dim * point_idx + k]
// << " != " << local_pos[k] << " " << point_idx);
// }
// }
// ++point_idx;
// }
// }
// DUMP("position coherency OK", DBG_INFO);
}
/* -------------------------------------------------------------------------- */
template <typename ContA, typename ContB>
void PointAssociation::updateForMigration(ContA &contA, ContB &contB) {
// constexpr UInt Dim = ContainerA::Dim;
DUMP("update migration for " << this->getID(), DBG_INFO);
DUMP("update local_points " << this->getID(), DBG_INFO);
STARTTIMER("syncMigration resize");
// update local number of points
if (this->in_group_A())
this->local_points = this->pointsAlist.size();
// update local number of points
if (this->in_group_B())
this->local_points = this->pointsBlist.size();
this->positions.resize(spatial_dimension * this->local_points);
STOPTIMER("syncMigration resize");
DUMP("synchMigration " << this->getID(), DBG_INFO);
STARTTIMER("syncMigration duo");
if ((this->in_group_A && contA.hasRefManager()) ||
(this->in_group_B && !contB.hasRefManager()))
this->synchronizeMigration(this->comm_group_A, this->comm_group_B);
if ((this->in_group_B && contB.hasRefManager()) ||
(this->in_group_A && !contA.hasRefManager()))
this->synchronizeMigration(this->comm_group_B, this->comm_group_A);
STOPTIMER("syncMigration duo");
DUMP("update release " << this->getID(), DBG_INFO);
LM_TOIMPLEMENT; // to review
// pointsAlist.setRelease(this->contA.getRelease());
// pointsBlist.setRelease(this->contB.getRelease());
if (check_coherency) {
DUMP("check position coherency " << this->getID(), DBG_INFO);
STARTTIMER("check position coherency");
this->checkPositionCoherency();
STOPTIMER("check position coherency");
}
}
/* -------------------------------------------------------------------------- */
/* LMDESC PointAssociation
This coupler computes the association between point degrees of freedom of
different domains
*/
/* LMHERITATE ComputeTrueDisplacement
*/
void PointAssociation::declareParams() {
DofAssociation::declareParams();
this->addSubParsableObject(true_disp);
/* LMKEYWORD TOL
max distance between two points to consider them coincident
*/
this->parseKeyword("TOL", tolerance);
/* LMKEYWORD PATH
file path to the filteratom file
*/
this->parseKeyword("PATH", path);
/* LMKEYWORD PBC
Specify if periodic boundary conditions should be considered for each
direction
*/
this->parseVectorKeyword("PBC", spatial_dimension, this->pbc);
/* LMKEYWORD CHECK_COHERENCY
Perform a systematic check of the communication scheme.
\textbf{Be careful, it is extremely computationally expensive}
*/
this->parseKeyword("CHECK_COHERENCY", check_coherency, false);
/* LMKEYWORD WEIGHT_PATH
file path to the python script to generate the weights
*/
this->parseKeyword("WEIGHT_PATH", weight_path);
}
/* -------------------------------------------------------------------------- */
// DECLARE_BRIDGING_ATOMIC_CONTINUUM_TEMPLATE(PointAssociation)
// DECLARE_BRIDGING_CONTINUUM_ATOMIC_TEMPLATE(PointAssociation)
__END_LIBMULTISCALE__
c4science · Help