Homec4science
Diffusion tamaas (features/netcdf-mpi) Browse

tamaas/src/contact_area.cpp
e96d4bebeabcfeatures/netcdf-mpi

Last Change
Raw File

contact_area.cpp
View Options

/**
*
* @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
*
* @section LICENSE
*
* Copyright (©) 2016 EPFL (Ecole Polytechnique Fédérale de
* Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des
* Solides)
*
* Tamaas 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.
*
* Tamaas 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 Tamaas. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#include"contact_area.hh"
#include <iostream>
#include "cluster_grow.hh"
/* -------------------------------------------------------------------------- */
__BEGIN_TAMAAS__
ContactArea::ContactArea(UInt n, Real L): Surface<int>(n,L) {
//num_cluster = 0;
}
/* -------------------------------------------------------------------------- */
template <>
ContactArea::ContactArea(const Map2dSquare<Real> & m) :
Surface<int>(m.size(),m.getL()){
UInt n = this->size();
for (UInt i = 0 ; i < n*n ; ++i)
this->at(i) = (m(i) > 0.);
}
/* -------------------------------------------------------------------------- */
ContactArea::~ContactArea() {}
/* -------------------------------------------------------------------------- */
void ContactArea::detectContactClusters() {
// starting state was push as 2 by vlad
ClusterGrow clg(*this,2);
UInt N = this->size();
for (UInt i = 0; i < N; ++i) {
for (UInt j = 0; j < N; ++j) {
clg.grow<true>(i,j);
}
}
for (UInt i = 0; i < N; ++i) {
for (UInt j = 0; j < N; ++j) {
if (this->at(i,j) == -1)
SURFACE_FATAL("problem " << i << " " << j << " is -1");
}
}
for (UInt i = 0; i < clg.getNbClusters(); ++i)
colors.push_back(i+3);
}
/* -------------------------------------------------------------------------- */
void ContactArea::replaceColors(int from, int to) {
// Repalce color "from" by the color "to"
UInt n = this->size();
for (UInt i = 0; i < n; i++)
for (UInt j = 0; j < n; j++)
if (this->at(i,j) == from)
this->at(i,j) = to;
}
/* -------------------------------------------------------------------------- */
void ContactArea::randomizeClusterID(UInt seed) {
std::cout << "/ Randomize colors of contact clusters." << std::endl;
UInt shift=1000, tryid;
bool found;
std::vector<UInt> cid;
UInt num_cluster = colors.size();
cid.resize(num_cluster);
while (shift < num_cluster)
shift*=10;
srand(seed);
for (int i = 0; i < int(colors.size()); i++) {
tryid=-1;
do {
tryid = shift+int((num_cluster+10)*random()/Real(RAND_MAX));
found = true;
for (int j = 0; j < i; j++)
if (cid[j] == tryid)
found = false;
}
while (!found);
cid[i]=tryid;
replaceColors(colors[i],tryid);
colors[i]=tryid;
}
// Shift all colors back to 1
for (int i = 0; i < int(colors.size()); i++) {
replaceColors(colors[i],colors[i]-shift+1);
colors[i]-=shift-1;
}
}
/* -------------------------------------------------------------------------- */
std::vector<ContactCluster> ContactArea::getContactClusters() {
std::vector<ContactCluster> clusters;
UInt nb_colors = this->colors.size();
for (UInt i = 0; i < nb_colors; i++) {
// std::cout << "Extracting contact cluster " << i << std::endl;
ContactCluster c(colors[i], *this);
//std::cout << "Compute statistics of contact cluster " << i << std::endl;
c.computeStatistics();
clusters.push_back(c);
}
return clusters;
}
/* -------------------------------------------------------------------------- */
__END_TAMAAS__
c4science · Help