Page MenuHomec4science
Files F65945473

bridging.cc
No OneTemporary
Actions

File Metadata

Created
Fri, Jun 7, 06:33

bridging.cc
View Options

/**
* @file bridging.cc
*
* @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
*
* @date Fri Jul 11 15:47:44 2014
*
* @brief Bridging object between atomistic and finite elements
*
* @section LICENSE
*
* Copyright INRIA and CEA
*
* The LibMultiScale is a C++ parallel framework for the multiscale
* coupling methods dedicated to material simulations. This framework
* provides an API which makes it possible to program coupled simulations
* and integration of already existing codes.
*
* This Project was initiated in a collaboration between INRIA Futurs Bordeaux
* within ScAlApplix team and CEA/DPTA Ile de France.
* The project is now continued at the Ecole Polytechnique Fédérale de Lausanne
* within the LSMS/ENAC laboratory.
*
* This software is governed by the CeCILL-C license under French law and
* abiding by the rules of distribution of free software. You can use,
* modify and/ or redistribute the software under the terms of the CeCILL-C
* license as circulated by CEA, CNRS and INRIA at the following URL
* "http://www.cecill.info".
*
* As a counterpart to the access to the source code and rights to copy,
* modify and redistribute granted by the license, users are provided only
* with a limited warranty and the software's author, the holder of the
* economic rights, and the successive licensors have only limited
* liability.
*
* In this respect, the user's attention is drawn to the risks associated
* with loading, using, modifying and/or developing or reproducing the
* software by the user in light of its specific status of free software,
* that may mean that it is complicated to manipulate, and that also
* therefore means that it is reserved for developers and experienced
* professionals having in-depth computer knowledge. Users are therefore
* encouraged to load and test the software's suitability as regards their
* requirements in conditions enabling the security of their systems and/or
* data to be ensured and, more generally, to use and operate it in the
* same conditions as regards security.
*
* The fact that you are presently reading this means that you have had
* knowledge of the CeCILL-C license and that you accept its terms.
*
*/
/* -------------------------------------------------------------------------- */
#include "bridging.hh"
#include "compute_extract.hh"
#include "factory_multiscale.hh"
#include "filter_geometry.hh"
#include "geometry_manager.hh"
#include "lm_common.hh"
#include "trace_atom.hh"
/* -------------------------------------------------------------------------- */
__BEGIN_LIBMULTISCALE__
Bridging::Bridging(const std::string &name)
: LMObject(name), DofAssociation(name),
unmatchedPointList("unmatched-point-" + name),
unmatchedMeshList("unmatched-fe-" + name),
pointList(createOutput("pointList")), meshList(createOutput("meshList")),
contMesh(createOutput("contMesh")), grid(createOutput("grid")) {
assoc_found = 0;
unmatchedPointList.setCommGroup(this->comm_group_A);
unmatchedMeshList.setCommGroup(this->comm_group_B);
// contains the elements selected by geometry
this->createOutput("unmatchedMeshList") = unmatchedMeshList;
// contains the points selected by geometry
this->createOutput("unmatchedPointList") = unmatchedPointList;
is_in_continuum = this->comm_group_B.amIinGroup();
is_in_atomic = this->comm_group_A.amIinGroup();
}
/* -------------------------------------------------------------------------- */
Bridging::~Bridging() {}
/* -------------------------------------------------------------------------- */
void Bridging::filterArray(ContainerArray<Real> &array) {
UInt index = 0;
for (UInt i = 0; i < pointToElement.size(); ++i)
if (pointToElement[i] != UINT_MAX) {
DUMP("for container, atom " << i << " becomes " << index, DBG_ALL);
LM_ASSERT(i < array.size() && index < array.size(),
"overflow detected: nmax = " << array.size() << " , index = "
<< index << " and i = " << i);
array(index) = array(i);
++index;
}
}
/* -------------------------------------------------------------------------- */
// void Bridging::cumulPBC(ContainerArray<Real> &data) {
// UInt npairs [[gnu::unused]] = local_pbc_pairs.size();
// DUMP("Detected " << npairs << " local pairs", DBG_INFO);
// for (auto &&pair : local_pbc_pairs) {
// UInt i1 = pair.first;
// UInt i2 = pair.second;
// data(i2) += data(i1);
// }
// for (auto &&pair : local_pbc_pairs) {
// UInt i1 = pair.first;
// UInt i2 = pair.second;
// data(i1) = data(i2);
// }
// }
/* -------------------------------------------------------------------------- */
void Bridging::mesh2Point(const ContainerArray<Real> &data_node,
ContainerArray<Real> &data_atom) {
if (this->local_points == 0) {
DUMP("Warning : atomic container is empty:"
<< " cannot interpolate mesh fields on control points (check "
"geometries ?!)",
DBG_WARNING);
return;
}
data_atom = smatrix.transpose() * data_node.matrix();
}
/* -------------------------------------------------------------------------- */
void Bridging::mesh2Point(FieldType field_type,
ContainerArray<Real> &data_point) {
ComputeExtract extracted_field("ComputeExtract:" + this->getID());
extracted_field.setParam("FIELD", field_type);
extracted_field.compute(this->meshList);
this->mesh2Point(extracted_field.evalArrayOutput(), data_point);
}
/* -------------------------------------------------------------------------- */
void Bridging::leastSquarePointData(ContainerArray<Real> &,
ContainerArray<Real> &, UInt) {
LM_TOIMPLEMENT;
}
/* -------------------------------------------------------------------------- */
void Bridging::copySlaveValues(ContainerArray<Real> &v) {
for (auto &&pair : local_pbc_pairs) {
UInt i1 = pair.first;
UInt i2 = pair.second;
v(i1) = v(i2);
}
}
/* -------------------------------------------------------------------------- */
/* Parallel Methods */
/* -------------------------------------------------------------------------- */
void Bridging::commPoint2Continuum(ContainerArray<Real> &buf) {
if (this->local_points != 0) {
this->distributeVectorA2B("communicate buffer", buf);
}
}
/* -------------------------------------------------------------------------- */
void Bridging::commContinuum2Point(ContainerArray<Real> &buf) {
if (this->local_points != 0) {
this->distributeVectorB2A("communicate buffer", buf);
}
}
/* -------------------------------------------------------------------------- */
void Bridging::synchSumBuffer(ContainerArray<Real> &buf) {
if (this->local_points != 0) {
this->synchronizeVectorBySum("synchsumvector", buf);
}
}
/* -------------------------------------------------------------------------- */
void Bridging::filterRejectedContinuumOwners(
std::vector<std::vector<UInt>> &unassociated_atoms) {
if (this->comm_group_A == this->comm_group_B)
return;
if (not this->in_group_B())
return;
DUMP("local points", DBG_MESSAGE);
if (this->local_points == 0)
return;
// loop through the unassociated atoms declared per each processor
// The total list of atom index rejected is constructed
UInt offset = 0;
UInt cpt_unassociated = 0;
std::vector<UInt> atom_to_remove;
for (auto &&[proc, c_w] : enumerate(com_with)) {
if (not c_w)
continue;
DUMP("atomic processor " << proc << " unassociated "
<< unassociated_atoms[proc].size() << " atoms ",
DBG_INFO_STARTUP);
for (auto &&unassociated : unassociated_atoms[proc]) {
// reconstruct the index for atom declared unassociated
UInt unassociated_index = offset + unassociated;
atom_to_remove.push_back(unassociated_index);
pointToElement[unassociated_index] = UINT_MAX;
// auto el_index = this->pointToElement[unassociated_index];
// LM_ASSERT(el_index != UINT_MAX,
// "this should not happend !! " << unassociated_index);
// // remove the column associated to atom in the shapematrix
// this->smatrix.removeAtom(unassociated_index);
// // alter the atom-element mapping
// pointToElement[unassociated_index] = UINT_MAX;
++cpt_unassociated;
DUMP("unassociated atom at position "
<< unassociated_index << " whereas offset = " << offset
<< " next offset is " << offset + this->nb_points_per_proc[proc],
DBG_INFO_STARTUP);
}
offset += this->nb_points_per_proc[proc];
}
smatrix.removeAtoms(atom_to_remove);
// corrects the local_points variable
if (cpt_unassociated) {
DUMP("removed " << cpt_unassociated
<< " atoms from unassociation. Now local atoms is "
<< this->local_points,
DBG_MESSAGE);
this->local_points -= cpt_unassociated;
DUMP("and becomes " << this->local_points, DBG_MESSAGE);
}
// compute contiguous indexes for atoms (compression of indexes)
// the result is a mapping from old indexes to new ones
// output: mapping_old_to_new_indexes
// std::vector<UInt> mapping_old_to_new_indexes;
// mapping_old_to_new_indexes.assign(pointToElement.size(), UINT_MAX);
// UInt indirection = 0;
// for (auto &&[at_index, el_index] : enumerate(pointToElement)) {
// if (el_index == UINT_MAX)
// continue;
// mapping_old_to_new_indexes[at_index] = indirection;
// ++indirection;
// }
// LM_ASSERT(indirection == this->local_points, "this should not happend "
// << indirection << " "
// << this->local_points);
// checks that number of mapped atoms equals local atoms variable
UInt cpt_tmp = 0;
for (auto &&el_index : pointToElement) {
if (el_index != UINT_MAX)
cpt_tmp++;
}
if (cpt_tmp != this->local_points)
LM_FATAL("biip this should not happend " << cpt_tmp << " "
<< this->local_points);
// // now renumber indirection in shapematrix
// indirection = 0;
// for (auto &&[at_index, el_index] : enumerate(pointToElement)) {
// if (el_index == UINT_MAX)
// continue;
// UInt new_index = mapping_old_to_new_indexes[at_index];
// LM_ASSERT(new_index != UINT_MAX, "This should not happend");
// smatrix.changeAtomIndirection(at_index, new_index, indirection);
// ++indirection;
// }
// this->smatrix.swapAtomIndirections();
// now i set the duo vector
this->createDuoVectorB("FE", pointToElement);
}
/* -------------------------------------------------------------------------- */
void Bridging::setPBCPairs(std::vector<std::pair<UInt, UInt>> &pairs) {
if (this->is_in_continuum)
pbc_pairs = pairs;
}
/* -------------------------------------------------------------------------- */
UInt Bridging::getNumberPoints() { LM_TOIMPLEMENT; }
/* -------------------------------------------------------------------------- */
UInt Bridging::getNumberElems() { LM_TOIMPLEMENT; }
/* -------------------------------------------------------------------------- */
UInt Bridging::getNumberNodes() { LM_TOIMPLEMENT; }
/* -------------------------------------------------------------------------- */
void Bridging::attachVector(ContainerArray<Real> &tab) {
this->pointList.get<ContainerInterface>().attachObject(tab);
}
/* -------------------------------------------------------------------------- */
void Bridging::updateForMigration() {
DUMPBYPROC("update migration for " << this->getID(), DBG_INFO, 0);
STARTTIMER("syncMigration resize");
// update local number of atoms in atomic part
if (this->in_group_A()) {
this->local_points =
this->getOutput<ContainerInterface>("pointList").size();
this->positions.resize(this->local_points, spatial_dimension);
}
STOPTIMER("syncMigration resize");
STARTTIMER("syncMigration duo");
this->synchronizeMigration(this->comm_group_A, this->comm_group_B);
STOPTIMER("syncMigration duo");
// pointList.setRelease(this->contA.getRelease());
// meshList.setRelease(this->contB.getRelease());
}
/* -------------------------------------------------------------------------- */
/* LMDESC Bridging
This class implements a bridging zone where
point and elements are associated.
For debugging pruposes, several set of DOf are registered to the central
system:
- unmatched-fe-\${couplerID} : the set of nodes and elements
that are contained in the geometry provided by GEOMETRY keyword.
- unmatched-point-\${couplerID} : the set of points
that are contained in the geometry provided by GEOMETRY keyword.
- matched-fe-\${couplerID} : the set of nodes and elements
that are containing at least one point of unmatched-point-\${couplerID}.
- matched-point-\${couplerID} : the set of points
that are contained in at least one element of unmatched-fe-\${couplerID}.
*/
/* LMHERITANCE filter_geometry dof_association */
void Bridging::declareParams() {
this->addSubParsableObject(unmatchedMeshList);
DofAssociation::declareParams();
}
/* -------------------------------------------------------------------------- */
void Bridging::setPointField(FieldType field, ContainerArray<Real> &buffer) {
StimulationField impose_field("impose_point_field");
impose_field.setParam("FIELD", field);
impose_field.compute("dofs"_input = this->pointList, "field"_input = buffer);
}
/* -------------------------------------------------------------------------- */
void Bridging::setPointField(FieldType field, Real val, UInt dim) {
ContainerArray<Real> data_point(pointList.get<ContainerInterface>().size(),
dim);
for (auto &&v : data_point.rowwise()) {
v = val;
}
setPointField(field, data_point);
}
/* -------------------------------------------------------------------------- */
void Bridging::addPointField(FieldType field, ContainerArray<Real> &buffer) {
StimulationField impose_field("impose_point_field");
impose_field.setParam("FIELD", field);
impose_field.setParam("ADDITIVE", true);
impose_field.compute("dofs"_input = this->pointList, "field"_input = buffer);
}
/* -------------------------------------------------------------------------- */
void Bridging::setMeshField(FieldType field, ContainerArray<Real> &buffer) {
StimulationField impose_field("impose_mesh_field");
impose_field.setParam("FIELD", field);
impose_field.compute("dofs"_input = this->meshList, "field"_input = buffer);
}
/* -------------------------------------------------------------------------- */
void Bridging::addMeshField(FieldType field, ContainerArray<Real> &buffer) {
StimulationField impose_field("impose_mesh_field");
impose_field.setParam("FIELD", field);
impose_field.setParam("ADDITIVE", true);
impose_field.compute("dofs"_input = this->meshList, "field"_input = buffer);
}
/* -------------------------------------------------------------------------- */
void Bridging::setPointField(FieldType field, Eigen::VectorXd val) {
ContainerArray<Real> data_point(pointList.get<ContainerInterface>().size(),
val.rows());
for (auto &&v : data_point.rowwise()) {
v = val;
}
setPointField(field, data_point);
}
/* -------------------------------------------------------------------------- */
void Bridging::setMeshField(FieldType field, Eigen::VectorXd val) {
ContainerArray<Real> data_mesh(pointList.get<ContainerInterface>().size(),
val.rows());
data_mesh = val;
setPointField(field, data_mesh);
}
/* -------------------------------------------------------------------------- */
__END_LIBMULTISCALE__

Event Timeline

c4science · Help