Homec4science
Diffusion LibMultiScale (features/new_parralelism) Browse

libmultiscale/src/md/lammps/import_lammps.cc
eabf556cc7ecfeatures/new_parralelism

Last Change
Raw File

import_lammps.cc
View Options

/**
* @file import_lammps.cc
*
* @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
* @author Jaehyun Cho <jaehyun.cho@epfl.ch>
*
* @date Wed Nov 20 07:29:13 2013
*
* @brief The programming hooks to let LM be aware of processor migrations
*
* @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 "lm_common.hh" #ifdef LIBMULTISCALE_USE_LAMMPS #include
// "domain_lammps.hh" #include "geometry.hh" #include
// "reference_manager_lammps.hh" #include "trace_atom.hh"
// /* --------------------------------------------------------------------------
// */ using namespace LAMMPS_NS; #include "atom_vec.h" #include
// "atom_vec_atomic.h" #include "domain.h" #include "fix.h" #include
// "lammps_common.hh" #include "memory.h"
// /* --------------------------------------------------------------------------
// */ #include "import_lammps.hh"
// /* --------------------------------------------------------------------------
// */ enum Topology {
// FACEXLEFT,
// FACEXRIGHT,
// FACEYLEFT,
// FACEYRIGHT,
// FACEZLEFT,
// FACEZRIGHT,
// EDGE1,
// EDGE2,
// EDGE3,
// EDGE4,
// EDGE5,
// EDGE6,
// EDGE7,
// EDGE8,
// EDGE9,
// EDGE10,
// EDGE11,
// EDGE12,
// CORNER1,
// CORNER2,
// CORNER3,
// CORNER4,
// CORNER5,
// CORNER6,
// CORNER7,
// CORNER8
// };
// /* --------------------------------------------------------------------------
// */ static const int n_neighs = 26;
// // static const int commPairs[n_neighs][2] = {
// // // faces pairs
// // {FACEXLEFT, FACEXRIGHT},
// // {FACEXRIGHT, FACEXLEFT},
// // {FACEYLEFT, FACEYRIGHT},
// // {FACEYRIGHT, FACEYLEFT},
// // {FACEZLEFT, FACEZRIGHT},
// // {FACEZRIGHT, FACEZLEFT},
// // // edges pairs
// // {EDGE1, EDGE2},
// // {EDGE2, EDGE1},
// // {EDGE3, EDGE4},
// // {EDGE4, EDGE3},
// // {EDGE5, EDGE6},
// // {EDGE6, EDGE5},
// // {EDGE7, EDGE8},
// // {EDGE8, EDGE7},
// // {EDGE9, EDGE10},
// // {EDGE10, EDGE9},
// // {EDGE11, EDGE12},
// // {EDGE12, EDGE11},
// // // corner pairs
// // {CORNER1, CORNER2},
// // {CORNER2, CORNER1},
// // {CORNER3, CORNER4},
// // {CORNER4, CORNER3},
// // {CORNER5, CORNER6},
// // {CORNER6, CORNER5},
// // {CORNER7, CORNER8},
// // {CORNER8, CORNER7}};
// /* --------------------------------------------------------------------------
// */ int topologyCoords[n_neighs][3] = {
// // faces
// {-1, 0, 0}, // FACEXLEFT
// {1, 0, 0}, // FACEXRIGHT
// {0, -1, 0}, // FACEYLEFT
// {0, 1, 0}, // FACEYRIGHT
// {0, 0, -1}, // FACEZLEFT
// {0, 0, 1}, // FACEZRIGHT
// // edges
// {1, 0, -1}, // EDGE1
// {-1, 0, -1}, // EDGE2
// {1, 0, 1}, // EDGE3
// {-1, 0, 1}, // EDGE4
// {-1, -1, 0}, // EDGE5
// {1, -1, 0}, // EDGE6
// {-1, 1, 0}, // EDGE7
// {1, 1, 0}, // EDGE8
// {0, -1, -1}, // EDGE9
// {0, 1, -1}, // EDGE1
// {0, -1, 1}, // EDGE1
// {0, 1, 1}, // EDGE1
// // corners
// {-1, -1, -1}, // CORNER1
// {-1, 1, -1}, // CORNER2
// {-1, 1, 1}, // CORNER3
// {-1, -1, 1}, // CORNER4
// {1, -1, -1}, // CORNER5
// {1, 1, -1}, // CORNER6
// {1, 1, 1}, // CORNER7
// {1, -1, 1} // CORNER8
// };
// const char topologyStrings[n_neighs][255] = {
// "FACEXLEFT", "FACEXRIGHT", "FACEYLEFT", "FACEYRIGHT", "FACEZLEFT",
// "FACEZRIGHT", "EDGE1", "EDGE2", "EDGE3", "EDGE4",
// "EDGE5", "EDGE6", "EDGE7", "EDGE8", "EDGE9",
// "EDGE10", "EDGE11", "EDGE12", "CORNER1", "CORNER2",
// "CORNER3", "CORNER4", "CORNER5", "CORNER6", "CORNER7",
// "CORNER8"};
// /* --------------------------------------------------------------------------
// */ #define stamp10 1000 #define stamp20 2000
// /* --------------------------------------------------------------------------
// */
// __BEGIN_LIBMULTISCALE__
// /* --------------------------------------------------------------------------
// */ void checkMPIComm(MPI_Status &status, int flag) {
// if (flag != MPI_SUCCESS || status.MPI_ERROR != MPI_SUCCESS) {
// char errStr[MPI_MAX_ERROR_STRING];
// int len;
// MPI_Error_string(status.MPI_ERROR, errStr, &len);
// LM_FATAL("MPI_ERROR = " << status.MPI_ERROR << "[" << errStr);
// }
// }
// /* --------------------------------------------------------------------------
// */
// template <UInt Dim>
// inline void ImportLammps<Dim>::init(LAMMPS_NS::Domain &domain
// [[gnu::unused]],
// LAMMPS_NS::Atom &atom [[gnu::unused]],
// LAMMPS_NS::Comm &comm [[gnu::unused]],
// MPI_Comm &world [[gnu::unused]],
// int me [[gnu::unused]],
// int nprocs [[gnu::unused]]) {
// // old_nlocal = atom.nlocal;
// // VIEW_ATOM(RefLammps<Dim>);
// // this->atom = &atom;
// // this->domain = &domain;
// // this->comm = &comm;
// // this->world = world;
// // this->me = me;
// // this->nprocs = nprocs;
// // setPeriodic(domain.xperiodic, domain.yperiodic, domain.zperiodic);
// // setupTopology();
// // MPI_Barrier(world);
// // DUMP("before exchange nlocal = " << atom.nlocal, DBG_INFO);
// // if (domain.triclinic == 0) {
// // Xlo = domain.sublo;
// // Xhi = domain.subhi;
// // } else {
// // Xlo = domain.sublo_lamda;
// // Xhi = domain.subhi_lamda;
// // }
// }
// /* --------------------------------------------------------------------------
// */
// // template <UInt Dim> inline void ImportLammps<Dim>::registerSends() {
// // // clear the buffers
// // for (UInt i = 0; i < buf_send.size(); ++i)
// // buf_send[i].clear();
// // int nlocal = atom->nlocal;
// // int proc[3];
// // for (UInt i = 0; i < 3; ++i)
// // DUMP("box dimension " << i << " " << Xlo[i] << " " << Xhi[i],
// DBG_INFO);
// // DUMP(" nlocal = " << nlocal, DBG_INFO);
// // int at = 0;
// // DUMP("doing loop over atoms for migration", DBG_INFO);
// // while (at < nlocal) {
// // Real *X = atom->x[at];
// // #ifndef LM_OPTIMIZED
// // Real *X0 = atom->x0[at];
// // Real *V = atom->v[at];
// // Real *F = atom->f[at];
// // #endif
// // // setting actual coordinates
// // for (UInt k = 0; k < 3; ++k) {
// // proc[k] = 0;
// // // search the new processor
// // if (X[k] < Xlo[k])
// // proc[k] = -1;
// // if (X[k] >= Xhi[k])
// // proc[k] = 1;
// // // special treatment for pbc
// // if (periodic[k]) {
// // if (domain->triclinic == 0) {
// // // test if i am at the lower border
// // if (Xlo[k] == domain->boxlo[k])
// // if (X[k] >= Xhi[k] && (X[k] - Xhi[k]) > (Xhi[k] - Xlo[k]))
// // proc[k] = -1;
// // // test if i am at the upper border
// // if (Xhi[k] == domain->boxhi[k])
// // if (X[k] < Xlo[k] && (Xlo[k] - X[k]) > (Xhi[k] - Xlo[k]))
// // proc[k] = 1;
// // } else {
// // // test if i am at the lower border
// // if (Xlo[k] == domain->boxlo_lamda[k])
// // if (X[k] >= Xhi[k] && (X[k] - Xhi[k]) > (Xhi[k] - Xlo[k]))
// // proc[k] = -1;
// // // test if i am at the upper border
// // if (Xhi[k] == domain->boxhi_lamda[k])
// // if (X[k] < Xlo[k] && (Xlo[k] - X[k]) > (Xhi[k] - Xlo[k]))
// // proc[k] = 1;
// // }
// // }
// // }
// // // DUMP("at: "<<at<<", nlocal: "<<nlocal<<", proc[0]: "<<proc[0]<<",
// // // proc[1]: "<<proc[1]<<", proc[2]: "<<proc[2],DBG_MESSAGE);
// // if (proc[0] == 0 && proc[1] == 0 && proc[2] == 0) {
// // DUMP("No proc change for this atom " << at, DBG_ALL);
// // ++at;
// // continue;
// // }
// // DUMP("migrating atom " << at << " nlocal << " << nlocal, DBG_ALL);
// // for (UInt k = 0; k < 3; ++k) {
// // DUMP("coord " << k << " " << X[k] << " " << X0[k] << " " << V[k] <<
// " "
// // << F[k] << " "
// // << "[" << Xlo[k] << " - " << Xhi[k] << "] " << proc[k]
// // << " " << proc[k],
// // DBG_ALL);
// // }
// // // DUMP("i: "<< at <<", X: " << X[0] << " is in proc " <<
// // // proc[0],DBG_MESSAGE);
// // // DUMP("i: "<< at <<", Y: " << X[1] << " is in proc " <<
// // // proc[1],DBG_MESSAGE);
// // // DUMP("i: "<< at <<", Z: " << X[2] << " is in proc " <<
// // // proc[2],DBG_MESSAGE);
// // // build the sending buffers
// // for (int i = 0; i < n_neighs; ++i) {
// // if (topologyCoords[i][0] == proc[0] && topologyCoords[i][1] ==
// proc[1]
// // &&
// // topologyCoords[i][2] == proc[2])
// // if (comProcs[i] != MPI_PROC_NULL)
// // registerSend(at, comProcs[i]);
// // }
// // // fill the hole made by leaving atom
// // atom->avec->copy(nlocal - 1, at, 0);
// // // moveAtom(nlocal-1,at);
// // // decrease total number of atoms
// // nlocal--;
// // }
// // atom->nlocal = nlocal;
// // DUMP("nlocal is now " << atom->nlocal, DBG_INFO);
// // }
// /* --------------------------------------------------------------------------
// */
// // template <UInt Dim>
// // inline void ImportLammps<Dim>::registerSend(UInt at, UInt toProc) {
// // #ifdef TRACE_ATOM
// // CHECK_TRACED(atom->x0[at]) {
// // VIEW_ATOM(RefLammps<Dim>);
// // DUMP(ATOM_OUT(atom->x0[at])
// // << " atom is leaving my proc for proc " << toProc,
// // DBG_MESSAGE);
// // DUMP(ATOM_OUT(atom->x0[at])
// // << " atom is leaving my proc for proc " << toProc,
// // DBG_WARNING);
// // internal_index = -1;
// // }
// // #endif
// // sendAtom(at, toProc);
// // UInt nsend = buf_send[toProc].size();
// // buf_send[toProc].resize(nsend + packSize);
// // DUMP("buf_send size is now " << buf_send[toProc].size(), DBG_DETAIL);
// // UInt added = atom->avec->pack_exchange(at, &buf_send[toProc][nsend]);
// // if (added != packSize)
// // LM_FATAL("internal error");
// //}
// /* --------------------------------------------------------------------------
// */
// template <UInt Dim> inline void ImportLammps<Dim>::unPackAtoms(UInt procDe) {
// UInt m = 0;
// std::vector<Real> &buf = buf_recv[procDe];
// DUMP("unpacking atom " << buf_recv[procDe].size() / packSize << " atoms",
// DBG_INFO);
// while (m < buf_recv[procDe].size()) {
// atom->avec->unpack_exchange(&buf[m]);
// acceptAtom(atom->nlocal - 1, procDe);
// {
// UInt i = atom->nlocal - 1;
// Real **ptr = atom->x;
// Real **x0ptr = atom->x0;
// Real **vptr = atom->v;
// SET_INTERNAL_TRACE_INDEX(x0ptr[i], i);
// if (CHECK_TRACED(x0ptr[i])) {
// DUMP("received traced atom from " << procDe, DBG_MESSAGE);
// VIEW_ATOM(RefLammps<Dim>);
// }
// if ((ptr[i][0] < Xlo[0] || ptr[i][0] >= Xhi[0]) || ptr[i][1] < Xlo[1]
// ||
// ptr[i][1] >= Xhi[1] || ptr[i][2] < Xlo[2] || ptr[i][2] >= Xhi[2]) {
// LM_FATAL("pb severe in migration process (received from "
// << procDe << " : atom i " << i << " x: " << ptr[i][0]
// << " xlo: " << Xlo[0] << " xhi: " << Xhi[0] << " "
// << " y: " << ptr[i][1] << " ylo: " << Xlo[1]
// << " yhi: " << Xhi[1] << " "
// << " z: " << ptr[i][2] << " zlo: " << Xlo[2]
// << " zhi: " << Xhi[2] << " "
// << " big trouble an atom went far far away cannot "
// << "control anything now : its position is " << x0ptr[i][0]
// << " " << x0ptr[i][1] << " " << x0ptr[i][2]
// << " and its velocity is " << vptr[i][0] << " " <<
// vptr[i][1]
// << " " << vptr[i][2]);
// }
// }
// }
// if (m != buf_recv[procDe].size())
// LM_FATAL("not all atoms where unpacked");
// buf_recv[procDe].clear();
// }
// /* --------------------------------------------------------------------------
// */
// template <UInt Dim> inline void ImportLammps<Dim>::checkConsistency() {
// MPI_Barrier(world);
// int test = 0;
// MPI_Allreduce(&(atom->nlocal), &test, 1, MPI_INT, MPI_SUM, world);
// DUMP("now nlocal = " << atom->nlocal << " nmax = " << atom->nmax,
// DBG_INFO);
// if (test != atom->natoms) {
// // ref_manager.printBilan();
// LM_FATAL("atoms were lost/created "
// << test << " != " << atom->natoms << " " << test - atom->natoms
// << std::endl
// << "# local atoms is " << atom->nlocal << std::endl
// << "old # local atoms is " << old_nlocal << std::endl
// << "abort");
// }
// if (atom->nmax < atom->nlocal)
// LM_FATAL("pb with allocation of the atoms");
// }
// /* --------------------------------------------------------------------------
// */
// // template <UInt Dim> inline void ImportLammps<Dim>::atomCommunications() {
// // //! pack the atoms needing a migration
// // // registerSends();
// // //! exchange the atoms
// // atomCommunication();
// // //! unpack the data
// // UInt ncom = proc_comm.size();
// // for (UInt i = 0; i < ncom; ++i)
// // unPackAtoms(proc_comm[i]);
// // // check constistency
// // checkConsistency();
// // }
// // /*
// --------------------------------------------------------------------------
// // */
// template <UInt Dim> inline void ImportLammps<Dim>::atomCommunication() {
// std::map<int, MPI_Request> requests;
// UInt ncom = proc_comm.size();
// //! asynchronous send phase
// for (UInt i = 0; i < ncom; ++i) {
// UInt toProc = proc_comm[i];
// if (int(toProc) == me)
// continue;
// int nbSend = buf_send[toProc].size();
// DUMP("Send to proc " << toProc << " the " << nbSend
// << " Reals t:" << stamp20 + me,
// DBG_INFO);
// // if (!nbSend) continue;
// if (requests.count(toProc))
// LM_FATAL("request to " << toProc << " already made");
// MPI_Request &req = requests[toProc];
// MPI_Isend(&buf_send[toProc][0], nbSend, MPI_DOUBLE, toProc, stamp20 + me,
// world, &req);
// }
// //! blocking receive phase
// for (UInt i = 0; i < ncom; ++i) {
// UInt fromProc = proc_comm[i];
// if (int(fromProc) == me)
// continue;
// MPI_Status status;
// status.MPI_ERROR = MPI_SUCCESS;
// DUMP("probing from proc " << fromProc, DBG_INFO);
// MPI_Probe(fromProc, stamp20 + fromProc, world, &status);
// int nbRecv = 0;
// MPI_Get_count(&status, MPI_DOUBLE, &nbRecv);
// buf_recv[fromProc].resize(nbRecv);
// DUMP("Receive from proc " << fromProc << " the " << nbRecv
// << " Reals t:" << stamp20 + fromProc,
// DBG_INFO);
// // if (!nbRecv) continue;
// int ret = MPI_Recv(&buf_recv[fromProc][0], nbRecv, MPI_DOUBLE, fromProc,
// stamp20 + fromProc, world, &status);
// checkMPIComm(status, ret);
// }
// //! wait for sends to be complete
// std::map<int, MPI_Request>::iterator it = requests.begin();
// std::map<int, MPI_Request>::iterator end = requests.end();
// while (it != end) {
// MPI_Status status;
// status.MPI_ERROR = MPI_SUCCESS;
// MPI_Request req = it->second;
// int toProc = it->first;
// DUMP("wait for com to " << toProc << " to be finished", DBG_INFO);
// int ret = MPI_Wait(&req, &status);
// if (buf_send[toProc].size()) {
// checkMPIComm(status, ret);
// }
// DUMP("send to proc " << toProc << " finished ", DBG_INFO);
// ++it;
// }
// }
// /* --------------------------------------------------------------------------
// */
// template <UInt Dim>
// int ImportLammps<Dim>::getProcNumber(int ind[3], int offset[3]) {
// /* Return the number of the processor located at the position
// i[0],i[1],i[2]
// */
// /* It returns MPI_PROC_NULL if none */
// /* It hanles periodic topology */
// int index[3];
// for (UInt i = 0; i < 3; ++i)
// index[i] = ind[i] + offset[i];
// DUMP("ind " << ind[0] << " " << ind[1] << " " << ind[2], DBG_DETAIL);
// DUMP("offset " << offset[0] << " " << offset[1] << " " << offset[2],
// DBG_DETAIL);
// DUMP("periodic " << periodic[0] << " " << periodic[1] << " " <<
// periodic[2],
// DBG_DETAIL);
// DUMP("nbrProc " << nbrProc[0] << " " << nbrProc[1] << " " << nbrProc[2],
// DBG_DETAIL);
// DUMP("index " << index[0] << " " << index[1] << " " << index[2],
// DBG_DETAIL);
// // exclude out of the domain faces (negative side)
// for (UInt i = 0; i < 3; ++i)
// if (index[i] < 0 && !periodic[i])
// return MPI_PROC_NULL;
// else if (index[i] < 0)
// index[i] = nbrProc[i] - 1;
// // exclude out of the domain faces (positive side)
// for (UInt i = 0; i < 3; ++i)
// if (index[i] >= nbrProc[i] && !periodic[i])
// return MPI_PROC_NULL;
// else if (index[i] >= nbrProc[i])
// index[i] = 0;
// DUMP("index " << index[0] << " " << index[1] << " " << index[2],
// DBG_DETAIL); return proc_matrix[index];
// }
// /* --------------------------------------------------------------------------
// */
// // //! setting processor topology with a stamp style
// // template <UInt Dim> void ImportLammps<Dim>::setupTopology() {
// // int reorder = 0;
// // // get number of processors in x,y and z directions
// // for (UInt i = 0; i < 3; ++i)
// // nbrProc[i] = comm->procgrid[i];
// // // MPI communicator in cartesian coordinates
// // MPI_Comm cartesian;
// // // creation of the topology
// // MPI_Cart_create(world, 3, nbrProc, periodic, reorder, &cartesian);
// // // fill the processor repartition matrix
// // proc_matrix.resize(nbrProc);
// // int coords[3];
// // for (int n = 0; n < nprocs; n++) {
// // MPI_Cart_coords(cartesian, n, 3, coords);
// // proc_matrix[coords] = n;
// // DUMP("proc " << n << " has cartesian coordinates : " << coords[0] << "
// "
// // << coords[1] << " " << coords[2],
// // DBG_INFO);
// // }
// // /* coordonnes du processeur */
// // MPI_Cart_coords(cartesian, me, 3, coords);
// // DUMP("my coordinates are " << coords[0] << " " << coords[1] << " "
// // << coords[2],
// // DBG_INFO);
// // MPI_Comm_free(&cartesian);
// // // initialize comprocs
// // for (int i = 0; i < n_neighs; ++i) {
// // DUMP("get " << topologyStrings[i] << " proc number", DBG_DETAIL);
// // comProcs[i] = getProcNumber(coords, topologyCoords[i]);
// // DUMP(topologyStrings[i] << " " << comProcs[i], DBG_INFO);
// // }
// // for (UInt i = 0; i < 3; ++i)
// // DUMP("nbrProc[" << i << "]" << nbrProc[i], DBG_INFO);
// // // initialize com buffers array
// // buf_send.resize(nprocs);
// // buf_recv.resize(nprocs);
// // // init list of processors to discuss with
// // std::set<int> proc_comm_with;
// // for (int i = 0; i < n_neighs; ++i)
// // proc_comm_with.insert(comProcs[i]);
// // proc_comm.clear();
// // std::set<int>::iterator it = proc_comm_with.begin();
// // std::set<int>::iterator end = proc_comm_with.end();
// // while (it != end) {
// // if (*it != MPI_PROC_NULL) {
// // proc_comm.push_back(*it);
// // DUMP("I should communicate with proc " << *it, DBG_DETAIL);
// // }
// // ++it;
// // }
// // DUMP("I should communicate with " << proc_comm.size() << " processors",
// // DBG_DETAIL);
// // // compute the packSize
// // double tmp[1000];
// // packSize = atom->avec->pack_exchange(0, tmp);
// // DUMP("packSize is " << packSize, DBG_INFO);
// // }
// /* --------------------------------------------------------------------------
// */ template class ImportLammps<2>; template class ImportLammps<3>;
// ImportLammpsInterface *ImportLammpsInterface::import = NULL;
// /* --------------------------------------------------------------------------
// */
// __END_LIBMULTISCALE__
// using namespace libmultiscale;
// /* ---------------------------------------------------------------------- */
// // void AtomVec::copy(int i, int j, int delflag) {
// // int m,n,datatype,cols,collength,ncols;
// // void *pdata,*plength;
// // tag[j] = tag[i];
// // type[j] = type[i];
// // mask[j] = mask[i];
// // image[j] = image[i];
// // x[j][0] = x[i][0];
// // x[j][1] = x[i][1];
// // x[j][2] = x[i][2];
// // x0[j][0] = x0[i][0];
// // x0[j][1] = x0[i][1];
// // x0[j][2] = x0[i][2];
// // v[j][0] = v[i][0];
// // v[j][1] = v[i][1];
// // v[j][2] = v[i][2];
// // if (ncopy) {
// // for (n = 0; n < ncopy; n++) {
// // pdata = mcopy.pdata[n];
// // datatype = mcopy.datatype[n];
// // cols = mcopy.cols[n];
// // if (datatype == Atom::DOUBLE) {
// // if (cols == 0) {
// // double *vec = *((double **) pdata);
// // vec[j] = vec[i];
// // } else if (cols > 0) {
// // double **array = *((double ***) pdata);
// // for (m = 0; m < cols; m++)
// // array[j][m] = array[i][m];
// // } else {
// // double **array = *((double ***) pdata);
// // collength = mcopy.collength[n];
// // plength = mcopy.plength[n];
// // if (collength) ncols = (*((int ***) plength))[i][collength-1];
// // else ncols = (*((int **) plength))[i];
// // for (m = 0; m < ncols; m++)
// // array[j][m] = array[i][m];
// // }
// // } else if (datatype == Atom::INT) {
// // if (cols == 0) {
// // int *vec = *((int **) pdata);
// // vec[j] = vec[i];
// // } else if (cols > 0) {
// // int **array = *((int ***) pdata);
// // for (m = 0; m < cols; m++)
// // array[j][m] = array[i][m];
// // } else {
// // int **array = *((int ***) pdata);
// // collength = mcopy.collength[n];
// // plength = mcopy.plength[n];
// // if (collength) ncols = (*((int ***) plength))[i][collength-1];
// // else ncols = (*((int **) plength))[i];
// // for (m = 0; m < ncols; m++)
// // array[j][m] = array[i][m];
// // }
// // } else if (datatype == Atom::BIGINT) {
// // if (cols == 0) {
// // bigint *vec = *((bigint **) pdata);
// // vec[j] = vec[i];
// // } else if (cols > 0) {
// // bigint **array = *((bigint ***) pdata);
// // for (m = 0; m < cols; m++)
// // array[j][m] = array[i][m];
// // } else {
// // bigint **array = *((bigint ***) pdata);
// // collength = mcopy.collength[n];
// // plength = mcopy.plength[n];
// // if (collength) ncols = (*((int ***) plength))[i][collength-1];
// // else ncols = (*((int **) plength))[i];
// // for (m = 0; m < ncols; m++)
// // array[j][m] = array[i][m];
// // }
// // }
// // }
// // }
// // if (bonus_flag) copy_bonus(i,j,delflag);
// // if (atom->nextra_grow)
// // for (int iextra = 0; iextra < atom->nextra_grow; iextra++)
// // modify->fix[atom->extra_grow[iextra]]->copy_arrays(i,j,delflag);
// // ImportLammpsInterface::getImport().moveAtom(i, j);
// // }
// /* --------------------------------------------------------------------------
// */
// // void Comm::exchange() {
// // MPI_Barrier(world);
// // DUMP("migration coms", DBG_INFO);
// // ImportLammpsInterface::getImport().init(*domain, *atom, *this, world,
// // this->me, this->nprocs);
// // // clear global->local map since atoms move & new ghosts are created
// // if (map_style)
// // atom->map_clear();
// // DUMP("after construction of sending buffers nlocal = "
// // << atom->nlocal << " nmax = " << atom->nmax,
// // DBG_INFO);
// // // send/receive communication to other procs
// // ImportLammpsInterface::getImport().atomCommunications();
// // MPI_Barrier(world);
// // ImportLammpsInterface::getImport().updateRefSubSets();
// // MPI_Barrier(world);
// // DUMP("migration all done for me !", DBG_INFO);
// // }
// /* --------------------------------------------------------------------------
// */
// // cleaned
// /* --------------------------------------------------------------------------
// */
// void AtomVec::create_atom(int itype, Real *coord) {
// // add the atom to my list of atoms
// if (!ImportLammpsInterface::getImport().getGeom(itype).contains(
// coord[0], coord[1], coord[2])) {
// return;
// }
// int m, n, datatype, cols;
// void *pdata;
// int nlocal = atom->nlocal;
// if (nlocal == nmax)
// grow(0);
// tag[nlocal] = 0;
// type[nlocal] = itype;
// x[nlocal][0] = coord[0];
// x[nlocal][1] = coord[1];
// x[nlocal][2] = coord[2];
// mask[nlocal] = 1;
// image[nlocal] =
// ((imageint)IMGMAX << IMG2BITS) | ((imageint)IMGMAX << IMGBITS) |
// IMGMAX;
// v[nlocal][0] = 0.0;
// v[nlocal][1] = 0.0;
// v[nlocal][2] = 0.0;
// // initialization additional fields
// for (n = 0; n < ncreate; n++) {
// pdata = mcreate.pdata[n];
// datatype = mcreate.datatype[n];
// cols = mcreate.cols[n];
// if (datatype == Atom::DOUBLE) {
// if (cols == 0) {
// double *vec = *((double **)pdata);
// vec[nlocal] = 0.0;
// } else {
// double **array = *((double ***)pdata);
// for (m = 0; m < cols; m++)
// array[nlocal][m] = 0.0;
// }
// } else if (datatype == Atom::INT) {
// if (cols == 0) {
// int *vec = *((int **)pdata);
// vec[nlocal] = 0;
// } else {
// int **array = *((int ***)pdata);
// for (m = 0; m < cols; m++)
// array[nlocal][m] = 0;
// }
// } else if (datatype == Atom::BIGINT) {
// if (cols == 0) {
// bigint *vec = *((bigint **)pdata);
// vec[nlocal] = 0;
// } else {
// bigint **array = *((bigint ***)pdata);
// for (m = 0; m < cols; m++)
// array[nlocal][m] = 0;
// }
// }
// }
// // if needed, initialize non-zero peratom values
// create_atom_post(nlocal);
// Real *X = x[nlocal];
// Real *xlo = domain->sublo;
// Real *xhi = domain->subhi;
// for (UInt i = 0; i < 3; ++i) {
// if (X[i] < xlo[i]) {
// LM_FATAL("cannot accept to create initially an atom outside"
// << " of the initial box"
// << "X[" << i << "] = " << X[i] << " < xlo = " << xlo[i]);
// }
// if (X[i] >= xhi[i]) {
// LM_FATAL("cannot accept to create initially an atom outside"
// << " of the initial box"
// << "X[" << i << "] = " << X[i] << " >= xhi = " << xhi[i]);
// }
// }
// atom->nlocal++;
// }
/* -------------------------------------------------------------------------- */
//#endif
c4science · Help