Homec4science
Diffusion LibMultiScale (features/pybind) Browse

libmultiscale/src/stimulation/stimulation_addforce.cc
b05d5f2aba78features/pybind

Last Change
Raw File

stimulation_addforce.cc
View Options

/**
* @file stimulation_addforce.cc
*
* @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
*
* @date Wed Jul 09 21:59:47 2014
*
* @brief Add a force to a set of atoms/nodes/points
*
* @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/>.
*
*/
/* -------------------------------------------------------------------------- */
#include "lm_common.hh"
#include "stimulation_addforce.hh"
#include "lib_md.hh"
#include "lib_dd.hh"
#include "lib_continuum.hh"
#include "filter.hh"
#include "ref_point_data.hh"
#include "communicator.hh"
/* -------------------------------------------------------------------------- */
__BEGIN_LIBMULTISCALE__
template <typename Cont>
void StimulationAddForce<Cont>::stimulate(Cont & cont, UInt stage) {
static const UInt Dim = Cont::Dim;
UInt cpt = cont.nbElem();
Communicator & comm = Communicator::getCommunicator();
CommGroup group = cont.getCommGroup();
UInt rank = comm.groupRank(lm_my_proc_id,group);
comm.allReduce(&cpt,1,group, "reduce total number of points concerned", OP_SUM);
if (rank == 0){
for (UInt i = 0; i < Dim; ++i)
DUMP("In direction " << i << " add an additional force "
<< total_force[i]/(Real)cpt + per_atom_force[i],DBG_INFO);
}
typename Cont::iterator it = cont.getIterator();
typedef typename Cont::Ref RefDof;
for (RefDof at = it.getFirst(); !it.end() ; at = it.getNext()) {
for (UInt i = 0; i < Dim; ++i) {
at.force(i) += total_force[i]/Real(cpt) + per_atom_force[i];
}
}
}
/* -------------------------------------------------------------------------- */
/* LMDESC ADDFORCE
Add a force to a set of atoms/nodes/points \\
The stimulator is provided a force $f^{req}$ and for each atoms/nodes/points in the
input the stimulator do:
$$ f += \frac{1}{\#points} f^{total} + f^{per-atom-force} $$
*/
/* LMEXAMPLE
STIMULATION force ADDFORCE TOTAL_FORCE 1 \\
STIMULATION force ADDFORCE PER_ATOM_FORCE 0.001
*/
template <typename Cont>
void StimulationAddForce<Cont>::declareParams(){
Stimulation<Cont>::declareParams();
this->changeDefault("STAGE", PRE_STEP3);
/* LMKEYWORD TOTAL_FORCE
Give the total force vector requested to be added. This quantity will
be divided by the number of atoms in the container.
*/
this->parseVectorKeyword("TOTAL_FORCE",spatial_dimension,total_force,VEC_DEFAULTS(0.,0.,0.));
/* LMKEYWORD PER_ATOM_FORCE
Give the force vector requested to be added on each atom.
*/
this->parseVectorKeyword("PER_ATOM_FORCE",spatial_dimension,per_atom_force,VEC_DEFAULTS(0.,0.,0.));
};
/* -------------------------------------------------------------------------- */
DECLARE_STIMULATION(StimulationAddForce,LIST_ATOM_MODEL)
DECLARE_STIMULATION(StimulationAddForce,LIST_CONTINUUM_MODEL)
DECLARE_STIMULATION(StimulationAddForce,LIST_DD_MODEL)
DECLARE_STIMULATION_REFPOINT(StimulationAddForce)
DECLARE_STIMULATION_GENERIC_MESH(StimulationAddForce)
__END_LIBMULTISCALE__
c4science · Help