Page MenuHomec4science
Files F67068432

dumper_lammps.cc
No OneTemporary
Actions

File Metadata

Created
Wed, Jun 19, 17:49

dumper_lammps.cc
View Options

/**
* @file dumper_lammps.cc
*
* @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
* @author Jaehyun Cho <jaehyun.cho@epfl.ch>
*
* @date Tue Aug 12 14:38:53 2014
*
* @brief This dumper allows to generate lammps text files
*
* @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 "dumper_lammps.hh"
#include "communicator.hh"
#include "compute_extract.hh"
#include "lib_continuum.hh"
#include "lib_dd.hh"
#include "lib_md.hh"
#include "lm_common.hh"
#include "ref_point_data.hh"
#include <fstream>
#include <iomanip>
/* -----------------------------------------------------------------*/
__BEGIN_LIBMULTISCALE__
template <typename Cont> void DumperLammps::dump(Cont &cont) {
constexpr UInt Dim = Cont::Dim;
// get the parallel context
CommGroup &group = cont.getCommGroup();
UInt my_rank = group.getMyRank();
UInt nb_procs = group.size();
// LM_ASSERT(nb_procs == 1, "DumperLammps is not yet parallel");
// prepare to retreive data from all other processors
ComputeExtract positions("Positions:ComputeExtract:" + this->getID());
positions.setParam("FIELD", _position);
positions.setCommGroup(cont.getCommGroup());
positions.build(cont);
ComputeExtract velocities("Velocity:ComputeExtract:" + this->getID());
velocities.setParam("FIELD", _velocity);
velocities.setCommGroup(cont.getCommGroup());
velocities.build(cont);
// ComputeExtract<Cont> tags("ComputeExtract:"+this->getID(),Dim);
// tags.setParam("FIELD",_tag);
// tags.build(cont);
// comm.reduce(&nb_data,1,group,"reduce total number of data element",OP_SUM);
// UInt nb_atoms = positions.getTotalNbData()/Dim;
UInt nb_atoms = positions.size();
group.reduce(&nb_atoms, 1, "reduce total number of data element", OP_SUM);
nb_atoms /= Dim;
std::ofstream file;
// open the file and write the headers
if (my_rank == 0) {
std::stringstream sstr;
sstr << this->getBaseName() << "_" << current_step << ".lammps";
file.open(sstr.str().c_str());
if (!file.is_open())
LM_FATAL("could not open file " << sstr.str());
file << "LAMMPS data file\n\n\n";
file << nb_atoms << " atoms\n\n";
UInt ntypes = geometries_for_types.size();
if (ntypes == 0)
++ntypes;
file << ntypes << " atom types\n\n";
const Vector<3> &xmin = cont.getBoundingBox().getXmin();
const Vector<3> &xmax = cont.getBoundingBox().getXmax();
for (UInt i = 0; i < 3; ++i) {
file << std::scientific << std::setprecision(15) << xmin[i] << " "
<< xmax[i] << " ";
if (i == 0)
file << "xlo xhi\n";
if (i == 1)
file << "ylo yhi\n";
if (i == 2)
file << "zlo zhi\n";
}
file << "Masses\n\n";
file << "1 " << std::scientific << std::setprecision(15)
<< cont.get(0).mass() << "\n\n";
// file << "Masses\n\n";
// file << "1 2.698154000000000e+01\n\n";
file << "Atoms\n\n";
}
UInt index_atom = 1;
for (UInt p = 0; p < nb_procs; ++p) {
std::vector<Real> data1;
UInt nb_data1 = positions.size();
if (my_rank == 0) {
if (p == 0) {
data1 = positions;
} else {
group.receive(&nb_data1, 1, p, "gatherData: receive number");
if (nb_data1 != 0) {
// data1.resize(nb_data);
group.receive(data1, p, "gatherData: receive data");
}
}
} else if (my_rank == p) {
group.send(&nb_data1, 1, 0, "gatherData: send number");
if (nb_data1 != 0) {
group.send(static_cast<std::vector<Real> &>(positions), 0,
"gatherData: send data");
}
}
// std::vector<Real> data3;
// UInt nb_data3 = tags.size();
// if (my_rank == 0){
// if(p == 0){
// data3 = tags;
// } else {
// comm.receive(&nb_data3,1,p,group,"gatherData: receive number");
// if (nb_data3 != 0) {
// //data1.resize(nb_data);
// comm.receive(data3,p,group,"gatherData: receive data");
// }
// }
// } else if (my_rank == p){
// comm.send(&nb_data3,1,0,group,"gatherData: send number");
// if (nb_data3 != 0){
// comm.send(static_cast<std::vector<Real>&
// >(tags),0,group,"gatherData: send data");
// }
// }
if (my_rank == 0) {
for (UInt j = 0; j < data1.size() / Dim; ++j, ++index_atom) {
Real X[3] = {0., 0., 0.};
for (UInt i = 0; i < Dim; ++i)
X[i] = data1[j * Dim + i];
UInt type = 1; // (UInt) data3[j];
for (UInt g = 0; g < geometries_for_types.size(); ++g) {
Geometry *geom = geometries_for_types[g];
if (geom->contains<Dim>(X))
type = g + 1;
}
file << index_atom << " " << type << " " << X[0] << " " << X[1] << " "
<< X[2] << std::endl;
}
}
}
if (my_rank == 0) {
file << "Velocities" << std::endl;
file << " " << std::endl;
}
index_atom = 1;
for (UInt p = 0; p < nb_procs; ++p) {
std::vector<Real> data2;
UInt nb_data = velocities.size();
if (my_rank == 0) {
if (p == 0) {
data2 = velocities;
} else {
group.receive(&nb_data, 1, p, "gatherData: receive number");
if (nb_data != 0) {
// data1.resize(nb_data);
group.receive(data2, p, "gatherData: receive data");
}
}
} else if (my_rank == p) {
group.send(&nb_data, 1, 0, "gatherData: send number");
if (nb_data != 0) {
group.send(static_cast<std::vector<Real> &>(velocities), 0,
"gatherData: send data");
}
}
if (my_rank == 0) {
for (UInt j = 0; j < data2.size() / Dim; ++j, ++index_atom) {
Real X[3] = {0., 0., 0.};
for (UInt i = 0; i < Dim; ++i)
X[i] = data2[j * Dim + i];
file << index_atom << " " << X[0] << " " << X[1] << " " << X[2]
<< std::endl;
}
}
}
file.close();
}
/* -------------------------------------------------------------------------- */
/* LMDESC LAMMPS
This dumper ouputs LAMMPS text files
*/
/* LMHERITANCE dumper */
/* LMEXAMPLE DUMPER lmp LAMMPS INPUT md FREQ 100 PREFIX ./ */
void DumperLammps::declareParams() {
DumperInterface::declareParams();
// Till commented this, since the variable is not used. In any case, the
// current position is dumped, never the initial
// this->parseKeyword("MODE",mode);
// parseCumulativeKeyword("ADDTYPE",ParserGeometry,geometries_for_types);
// if ("ADDTYPE" == key){
// GeomID geom;
// Parser::parseGeometry(geom,line);
// geometries_for_types.push_back(GeometryManager::getManager().getGeometry(geom));
// return true;
// }
}
/* -------------------------------------------------------------------------- */
DumperLammps::DumperLammps(const std::string &name) : LMObject(name) {
mode = 0;
}
/* -------------------------------------------------------------------------- */
DumperLammps::~DumperLammps() {}
/* -------------------------------------------------------------------------- */
DECLARE_DUMPER_MAKE_CALL(DumperLammps);
__END_LIBMULTISCALE__

Event Timeline

c4science · Help