Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-2018 The plumed team
3 : (see the PEOPLE file at the root of the distribution for a list of names)
4 :
5 : See http://www.plumed.org for more information.
6 :
7 : This file is part of plumed, version 2.
8 :
9 : plumed is free software: you can redistribute it and/or modify
10 : it under the terms of the GNU Lesser General Public License as published by
11 : the Free Software Foundation, either version 3 of the License, or
12 : (at your option) any later version.
13 :
14 : plumed is distributed in the hope that it will be useful,
15 : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 : GNU Lesser General Public License for more details.
18 :
19 : You should have received a copy of the GNU Lesser General Public License
20 : along with plumed. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 : #include "Atoms.h"
23 : #include "ActionAtomistic.h"
24 : #include "MDAtoms.h"
25 : #include "PlumedMain.h"
26 : #include "tools/Pbc.h"
27 : #include <algorithm>
28 : #include <iostream>
29 : #include <string>
30 : #include <cmath>
31 :
32 : using namespace std;
33 :
34 : namespace PLMD {
35 :
36 : /// We assume that charges and masses are constant along the simulation
37 : /// Set this to false if you want to revert to the original (expensive) behavior
38 : static const bool shareMassAndChargeOnlyAtFirstStep=true;
39 :
40 : class PlumedMain;
41 :
42 1159 : Atoms::Atoms(PlumedMain&plumed):
43 : natoms(0),
44 1159 : pbc(*new Pbc),
45 : md_energy(0.0),
46 : energy(0.0),
47 : dataCanBeSet(false),
48 : collectEnergy(false),
49 : energyHasBeenSet(false),
50 : positionsHaveBeenSet(0),
51 : massesHaveBeenSet(false),
52 : chargesHaveBeenSet(false),
53 : boxHasBeenSet(false),
54 : forcesHaveBeenSet(0),
55 : virialHasBeenSet(false),
56 : massAndChargeOK(false),
57 : shuffledAtoms(0),
58 : plumed(plumed),
59 : naturalUnits(false),
60 : MDnaturalUnits(false),
61 : timestep(0.0),
62 : forceOnEnergy(0.0),
63 : zeroallforces(false),
64 : kbT(0.0),
65 : asyncSent(false),
66 : atomsNeeded(false),
67 2318 : ddStep(0)
68 : {
69 1159 : mdatoms=MDAtomsBase::create(sizeof(double));
70 1159 : }
71 :
72 2318 : Atoms::~Atoms() {
73 1159 : if(actions.size()>0) {
74 0 : std::cerr<<"WARNING: there is some inconsistency in action added to atoms, as some of them were not properly destroyed. This might indicate an internal bug!!\n";
75 : }
76 1159 : delete mdatoms;
77 1159 : delete &pbc;
78 1159 : }
79 :
80 21537 : void Atoms::startStep() {
81 21537 : collectEnergy=false; energyHasBeenSet=false; positionsHaveBeenSet=0;
82 21537 : massesHaveBeenSet=false; chargesHaveBeenSet=false; boxHasBeenSet=false;
83 21537 : forcesHaveBeenSet=0; virialHasBeenSet=false; dataCanBeSet=true;
84 21537 : }
85 :
86 20963 : void Atoms::setBox(void*p) {
87 20963 : mdatoms->setBox(p);
88 20963 : Tensor b; mdatoms->getBox(b); boxHasBeenSet=true;
89 20963 : }
90 :
91 20991 : void Atoms::setPositions(void*p) {
92 20991 : plumed_massert( dataCanBeSet,"setPositions must be called after setStep in MD code interface");
93 20991 : plumed_massert( p || gatindex.size()==0, "NULL position pointer with non-zero local atoms");
94 20991 : mdatoms->setp(p); positionsHaveBeenSet=3;
95 20991 : }
96 :
97 20991 : void Atoms::setMasses(void*p) {
98 20991 : plumed_massert( dataCanBeSet,"setMasses must be called after setStep in MD code interface");
99 20991 : plumed_massert( p || gatindex.size()==0, "NULL mass pointer with non-zero local atoms");
100 20991 : mdatoms->setm(p); massesHaveBeenSet=true;
101 :
102 20991 : }
103 :
104 18773 : void Atoms::setCharges(void*p) {
105 18773 : plumed_massert( dataCanBeSet, "setCharges must be called after setStep in MD code interface");
106 18773 : plumed_massert( p || gatindex.size()==0, "NULL charges pointer with non-zero local atoms");
107 18773 : mdatoms->setc(p); chargesHaveBeenSet=true;
108 18773 : }
109 :
110 18841 : void Atoms::setVirial(void*p) {
111 18841 : plumed_massert( dataCanBeSet,"setVirial must be called after setStep in MD code interface");
112 18841 : mdatoms->setVirial(p); virialHasBeenSet=true;
113 :
114 18841 : }
115 :
116 2150 : void Atoms::setEnergy(void*p) {
117 2150 : plumed_massert( dataCanBeSet,"setEnergy must be called after setStep in MD code interface");
118 2150 : MD2double(p,md_energy);
119 2150 : md_energy*=MDUnits.getEnergy()/units.getEnergy();
120 2150 : energyHasBeenSet=true;
121 2150 : }
122 :
123 20991 : void Atoms::setForces(void*p) {
124 20991 : plumed_massert( dataCanBeSet,"setForces must be called after setStep in MD code interface");
125 20991 : plumed_massert( p || gatindex.size()==0, "NULL force pointer with non-zero local atoms");
126 20991 : forcesHaveBeenSet=3;
127 20991 : mdatoms->setf(p);
128 20991 : }
129 :
130 0 : void Atoms::setPositions(void*p,int i) {
131 0 : plumed_massert( dataCanBeSet,"setPositions must be called after setStep in MD code interface");
132 0 : plumed_massert( p || gatindex.size()==0, "NULL positions pointer with non-zero local atoms");
133 0 : mdatoms->setp(p,i); positionsHaveBeenSet++;
134 0 : }
135 :
136 0 : void Atoms::setForces(void*p,int i) {
137 0 : plumed_massert( dataCanBeSet,"setForces must be called after setStep in MD code interface");
138 0 : plumed_massert( p || gatindex.size()==0, "NULL force pointer with non-zero local atoms");
139 0 : mdatoms->setf(p,i); forcesHaveBeenSet++;
140 0 : }
141 :
142 19518 : void Atoms::share() {
143 19518 : std::set<AtomNumber> unique;
144 : // At first step I scatter all the atoms so as to store their mass and charge
145 : // Notice that this works with the assumption that charges and masses are
146 : // not changing during the simulation!
147 19518 : if(!massAndChargeOK && shareMassAndChargeOnlyAtFirstStep) {
148 322 : shareAll();
149 19840 : return;
150 : }
151 79111 : for(unsigned i=0; i<actions.size(); i++) if(actions[i]->isActive()) {
152 54392 : if(dd && shuffledAtoms>0) {
153 1132 : unique.insert(actions[i]->getUnique().begin(),actions[i]->getUnique().end());
154 : }
155 54392 : if(!actions[i]->getUnique().empty()) atomsNeeded=true;
156 : }
157 19196 : share(unique);
158 : }
159 :
160 394 : void Atoms::shareAll() {
161 394 : std::set<AtomNumber> unique;
162 394 : if(dd && shuffledAtoms>0)
163 106 : for(int i=0; i<natoms; i++) unique.insert(AtomNumber::index(i));
164 394 : atomsNeeded=true;
165 394 : share(unique);
166 394 : }
167 :
168 19590 : void Atoms::share(const std::set<AtomNumber>& unique) {
169 19590 : plumed_assert( positionsHaveBeenSet==3 && massesHaveBeenSet );
170 19590 : virial.zero();
171 19590 : if(zeroallforces || int(gatindex.size())==natoms) {
172 19314 : for(int i=0; i<natoms; i++) forces[i].zero();
173 : } else {
174 276 : for(unsigned i=0; i<gatindex.size(); i++) forces[gatindex[i]].zero();
175 : }
176 19590 : for(unsigned i=getNatoms(); i<positions.size(); i++) forces[i].zero(); // virtual atoms
177 19590 : forceOnEnergy=0.0;
178 19590 : mdatoms->getBox(box);
179 :
180 39180 : if(!atomsNeeded) return;
181 :
182 18883 : atomsNeeded=false;
183 :
184 18883 : if(int(gatindex.size())==natoms && shuffledAtoms==0) {
185 : // faster version, which retrieves all atoms
186 18559 : mdatoms->getPositions(0,natoms,positions);
187 : } else {
188 : // version that picks only atoms that are available on this proc
189 324 : mdatoms->getPositions(gatindex,positions);
190 : }
191 : // how many double per atom should be scattered:
192 18883 : int ndata=3;
193 18883 : if(!massAndChargeOK) {
194 322 : ndata=5;
195 322 : masses.assign(masses.size(),NAN);
196 322 : charges.assign(charges.size(),NAN);
197 322 : mdatoms->getCharges(gatindex,charges);
198 322 : mdatoms->getMasses(gatindex,masses);
199 : }
200 :
201 18883 : if(dd && shuffledAtoms>0) {
202 324 : if(dd.async) {
203 304 : for(unsigned i=0; i<dd.mpi_request_positions.size(); i++) dd.mpi_request_positions[i].wait();
204 304 : for(unsigned i=0; i<dd.mpi_request_index.size(); i++) dd.mpi_request_index[i].wait();
205 : }
206 324 : int count=0;
207 26246 : for(std::set<AtomNumber>::const_iterator p=unique.begin(); p!=unique.end(); ++p) {
208 25922 : if(dd.g2l[p->index()]>=0) {
209 11765 : dd.indexToBeSent[count]=p->index();
210 11765 : dd.positionsToBeSent[ndata*count+0]=positions[p->index()][0];
211 11765 : dd.positionsToBeSent[ndata*count+1]=positions[p->index()][1];
212 11765 : dd.positionsToBeSent[ndata*count+2]=positions[p->index()][2];
213 11765 : if(!massAndChargeOK) {
214 1245 : dd.positionsToBeSent[ndata*count+3]=masses[p->index()];
215 1245 : dd.positionsToBeSent[ndata*count+4]=charges[p->index()];
216 : }
217 11765 : count++;
218 : }
219 : }
220 324 : if(dd.async) {
221 304 : asyncSent=true;
222 304 : dd.mpi_request_positions.resize(dd.Get_size());
223 304 : dd.mpi_request_index.resize(dd.Get_size());
224 992 : for(int i=0; i<dd.Get_size(); i++) {
225 688 : dd.mpi_request_index[i]=dd.Isend(&dd.indexToBeSent[0],count,i,666);
226 688 : dd.mpi_request_positions[i]=dd.Isend(&dd.positionsToBeSent[0],ndata*count,i,667);
227 : }
228 : } else {
229 20 : const int n=(dd.Get_size());
230 20 : vector<int> counts(n);
231 40 : vector<int> displ(n);
232 40 : vector<int> counts5(n);
233 40 : vector<int> displ5(n);
234 20 : dd.Allgather(count,counts);
235 20 : displ[0]=0;
236 20 : for(int i=1; i<n; ++i) displ[i]=displ[i-1]+counts[i-1];
237 20 : for(int i=0; i<n; ++i) counts5[i]=counts[i]*ndata;
238 20 : for(int i=0; i<n; ++i) displ5[i]=displ[i]*ndata;
239 20 : dd.Allgatherv(&dd.indexToBeSent[0],count,&dd.indexToBeReceived[0],&counts[0],&displ[0]);
240 20 : dd.Allgatherv(&dd.positionsToBeSent[0],ndata*count,&dd.positionsToBeReceived[0],&counts5[0],&displ5[0]);
241 20 : int tot=displ[n-1]+counts[n-1];
242 1620 : for(int i=0; i<tot; i++) {
243 1600 : positions[dd.indexToBeReceived[i]][0]=dd.positionsToBeReceived[ndata*i+0];
244 1600 : positions[dd.indexToBeReceived[i]][1]=dd.positionsToBeReceived[ndata*i+1];
245 1600 : positions[dd.indexToBeReceived[i]][2]=dd.positionsToBeReceived[ndata*i+2];
246 1600 : if(!massAndChargeOK) {
247 432 : masses[dd.indexToBeReceived[i]] =dd.positionsToBeReceived[ndata*i+3];
248 432 : charges[dd.indexToBeReceived[i]] =dd.positionsToBeReceived[ndata*i+4];
249 : }
250 20 : }
251 : }
252 : }
253 : }
254 :
255 19590 : void Atoms::wait() {
256 19590 : dataCanBeSet=false; // Everything should be set by this stage
257 : // How many double per atom should be scattered
258 19590 : int ndata=3;
259 19590 : if(!massAndChargeOK)ndata=5;
260 :
261 19590 : if(dd) {
262 8274 : dd.Bcast(box,0);
263 : }
264 19590 : pbc.setBox(box);
265 :
266 19590 : if(collectEnergy) energy=md_energy;
267 :
268 19590 : if(dd && shuffledAtoms>0) {
269 : // receive toBeReceived
270 324 : if(asyncSent) {
271 : Communicator::Status status;
272 304 : int count=0;
273 992 : for(int i=0; i<dd.Get_size(); i++) {
274 688 : dd.Recv(&dd.indexToBeReceived[count],dd.indexToBeReceived.size()-count,i,666,status);
275 688 : int c=status.Get_count<int>();
276 688 : dd.Recv(&dd.positionsToBeReceived[ndata*count],dd.positionsToBeReceived.size()-ndata*count,i,667);
277 688 : count+=c;
278 : }
279 24626 : for(int i=0; i<count; i++) {
280 24322 : positions[dd.indexToBeReceived[i]][0]=dd.positionsToBeReceived[ndata*i+0];
281 24322 : positions[dd.indexToBeReceived[i]][1]=dd.positionsToBeReceived[ndata*i+1];
282 24322 : positions[dd.indexToBeReceived[i]][2]=dd.positionsToBeReceived[ndata*i+2];
283 24322 : if(!massAndChargeOK) {
284 2706 : masses[dd.indexToBeReceived[i]] =dd.positionsToBeReceived[ndata*i+3];
285 2706 : charges[dd.indexToBeReceived[i]] =dd.positionsToBeReceived[ndata*i+4];
286 : }
287 : }
288 304 : asyncSent=false;
289 : }
290 324 : if(collectEnergy) dd.Sum(energy);
291 : }
292 : // I take note that masses and charges have been set once for all
293 : // at the beginning of the simulation.
294 19590 : if(shareMassAndChargeOnlyAtFirstStep) massAndChargeOK=true;
295 19590 : }
296 :
297 19518 : void Atoms::updateForces() {
298 19518 : plumed_assert( forcesHaveBeenSet==3 );
299 19518 : if(forceOnEnergy*forceOnEnergy>epsilon) {
300 50 : double alpha=1.0-forceOnEnergy;
301 50 : mdatoms->rescaleForces(gatindex,alpha);
302 : }
303 19518 : mdatoms->updateForces(gatindex,forces);
304 19518 : if( !plumed.novirial && dd.Get_rank()==0 ) {
305 14811 : plumed_assert( virialHasBeenSet );
306 14811 : mdatoms->updateVirial(virial);
307 : }
308 19518 : }
309 :
310 331 : void Atoms::setNatoms(int n) {
311 331 : natoms=n;
312 331 : positions.resize(n);
313 331 : forces.resize(n);
314 331 : masses.resize(n);
315 331 : charges.resize(n);
316 331 : gatindex.resize(n);
317 331 : for(unsigned i=0; i<gatindex.size(); i++) gatindex[i]=i;
318 331 : }
319 :
320 :
321 1356 : void Atoms::add(const ActionAtomistic*a) {
322 1356 : actions.push_back(a);
323 1356 : }
324 :
325 1356 : void Atoms::remove(const ActionAtomistic*a) {
326 1356 : vector<const ActionAtomistic*>::iterator f=find(actions.begin(),actions.end(),a);
327 1356 : plumed_massert(f!=actions.end(),"cannot remove an action registered to atoms");
328 1356 : actions.erase(f);
329 1356 : }
330 :
331 :
332 115 : void Atoms::DomainDecomposition::enable(Communicator& c) {
333 115 : on=true;
334 115 : Set_comm(c.Get_comm());
335 115 : async=Get_size()<10;
336 115 : if(std::getenv("PLUMED_ASYNC_SHARE")) {
337 4 : std::string s(std::getenv("PLUMED_ASYNC_SHARE"));
338 4 : if(s=="yes") async=true;
339 4 : else if(s=="no") async=false;
340 0 : else plumed_merror("PLUMED_ASYNC_SHARE variable is set to " + s + "; should be yes or no");
341 : }
342 115 : }
343 :
344 265 : void Atoms::setAtomsNlocal(int n) {
345 265 : gatindex.resize(n);
346 265 : if(dd) {
347 265 : dd.g2l.resize(natoms,-1);
348 : // Since these vectors are sent with MPI by using e.g.
349 : // &dd.positionsToBeSent[0]
350 : // we make sure they are non-zero-sized so as to
351 : // avoid errors when doing boundary check
352 265 : if(n==0) n++;
353 265 : dd.positionsToBeSent.resize(n*5,0.0);
354 265 : dd.positionsToBeReceived.resize(natoms*5,0.0);
355 265 : dd.indexToBeSent.resize(n,0);
356 265 : dd.indexToBeReceived.resize(natoms,0);
357 : };
358 265 : }
359 :
360 134 : void Atoms::setAtomsGatindex(int*g,bool fortran) {
361 134 : plumed_massert( g || gatindex.size()==0, "NULL gatindex pointer with non-zero local atoms");
362 134 : ddStep=plumed.getStep();
363 134 : if(fortran) {
364 0 : for(unsigned i=0; i<gatindex.size(); i++) gatindex[i]=g[i]-1;
365 : } else {
366 134 : for(unsigned i=0; i<gatindex.size(); i++) gatindex[i]=g[i];
367 : }
368 134 : for(unsigned i=0; i<dd.g2l.size(); i++) dd.g2l[i]=-1;
369 134 : if( gatindex.size()==natoms ) {
370 0 : shuffledAtoms=0;
371 0 : for(unsigned i=0; i<gatindex.size(); i++) {
372 0 : if( gatindex[i]!=i ) { shuffledAtoms=1; break; }
373 : }
374 : } else {
375 134 : shuffledAtoms=1;
376 : }
377 134 : if(dd) {
378 134 : dd.Sum(shuffledAtoms);
379 134 : for(unsigned i=0; i<gatindex.size(); i++) dd.g2l[gatindex[i]]=i;
380 : }
381 134 : }
382 :
383 131 : void Atoms::setAtomsContiguous(int start) {
384 131 : ddStep=plumed.getStep();
385 131 : for(unsigned i=0; i<gatindex.size(); i++) gatindex[i]=start+i;
386 131 : for(unsigned i=0; i<dd.g2l.size(); i++) dd.g2l[i]=-1;
387 131 : if(dd) for(unsigned i=0; i<gatindex.size(); i++) dd.g2l[gatindex[i]]=i;
388 131 : if(gatindex.size()<natoms) shuffledAtoms=1;
389 131 : }
390 :
391 292 : void Atoms::setRealPrecision(int p) {
392 292 : delete mdatoms;
393 292 : mdatoms=MDAtomsBase::create(p);
394 292 : }
395 :
396 331 : int Atoms::getRealPrecision()const {
397 331 : return mdatoms->getRealPrecision();
398 : }
399 :
400 3310 : void Atoms::MD2double(const void*m,double&d)const {
401 3310 : plumed_assert(mdatoms); mdatoms->MD2double(m,d);
402 3310 : }
403 191 : void Atoms::double2MD(const double&d,void*m)const {
404 191 : plumed_assert(mdatoms); mdatoms->double2MD(d,m);
405 191 : }
406 :
407 331 : void Atoms::updateUnits() {
408 331 : mdatoms->setUnits(units,MDUnits);
409 331 : }
410 :
411 292 : void Atoms::setTimeStep(void*p) {
412 292 : MD2double(p,timestep);
413 292 : }
414 :
415 676533 : double Atoms::getTimeStep()const {
416 676533 : return timestep/units.getTime()*MDUnits.getTime();
417 : }
418 :
419 4 : void Atoms::setKbT(void*p) {
420 4 : MD2double(p,kbT);
421 4 : }
422 :
423 417 : double Atoms::getKbT()const {
424 417 : return kbT/units.getEnergy()*MDUnits.getEnergy();
425 : }
426 :
427 :
428 5 : void Atoms::createFullList(int*n) {
429 5 : vector<AtomNumber> fullListTmp;
430 15 : for(unsigned i=0; i<actions.size(); i++) if(actions[i]->isActive())
431 5 : fullListTmp.insert(fullListTmp.end(),actions[i]->getUnique().begin(),actions[i]->getUnique().end());
432 5 : std::sort(fullListTmp.begin(),fullListTmp.end());
433 5 : int nn=std::unique(fullListTmp.begin(),fullListTmp.end())-fullListTmp.begin();
434 5 : fullList.resize(nn);
435 5 : for(int i=0; i<nn; ++i) fullList[i]=fullListTmp[i].index();
436 5 : *n=nn;
437 5 : }
438 :
439 5 : void Atoms::getFullList(int**x) {
440 5 : if(!fullList.empty()) *x=&fullList[0];
441 1 : else *x=NULL;
442 5 : }
443 :
444 5 : void Atoms::clearFullList() {
445 5 : fullList.resize(0);
446 5 : }
447 :
448 331 : void Atoms::init() {
449 : // Default: set domain decomposition to NO-decomposition, waiting for
450 : // further instruction
451 331 : if(dd) {
452 115 : setAtomsNlocal(natoms);
453 115 : setAtomsContiguous(0);
454 : }
455 331 : }
456 :
457 115 : void Atoms::setDomainDecomposition(Communicator& comm) {
458 115 : dd.enable(comm);
459 115 : }
460 :
461 234 : void Atoms::resizeVectors(unsigned n) {
462 234 : positions.resize(n);
463 234 : forces.resize(n);
464 234 : masses.resize(n);
465 234 : charges.resize(n);
466 234 : }
467 :
468 117 : AtomNumber Atoms::addVirtualAtom(ActionWithVirtualAtom*a) {
469 117 : unsigned n=positions.size();
470 117 : resizeVectors(n+1);
471 117 : virtualAtomsActions.push_back(a);
472 117 : return AtomNumber::index(n);
473 : }
474 :
475 117 : void Atoms::removeVirtualAtom(ActionWithVirtualAtom*a) {
476 117 : unsigned n=positions.size();
477 117 : plumed_massert(a==virtualAtomsActions[virtualAtomsActions.size()-1],"virtual atoms should be destroyed in reverse creation order");
478 117 : resizeVectors(n-1);
479 117 : virtualAtomsActions.pop_back();
480 117 : }
481 :
482 63 : void Atoms::insertGroup(const std::string&name,const std::vector<AtomNumber>&a) {
483 63 : plumed_massert(groups.count(name)==0,"group named "+name+" already exists");
484 63 : groups[name]=a;
485 63 : }
486 :
487 63 : void Atoms::removeGroup(const std::string&name) {
488 63 : plumed_massert(groups.count(name)==1,"cannot remove group named "+name);
489 63 : groups.erase(name);
490 63 : }
491 :
492 72 : void Atoms::writeBinary(std::ostream&o)const {
493 72 : o.write(reinterpret_cast<const char*>(&positions[0][0]),natoms*3*sizeof(double));
494 72 : o.write(reinterpret_cast<const char*>(&box(0,0)),9*sizeof(double));
495 72 : o.write(reinterpret_cast<const char*>(&energy),sizeof(double));
496 72 : }
497 :
498 72 : void Atoms::readBinary(std::istream&i) {
499 72 : i.read(reinterpret_cast<char*>(&positions[0][0]),natoms*3*sizeof(double));
500 72 : i.read(reinterpret_cast<char*>(&box(0,0)),9*sizeof(double));
501 72 : i.read(reinterpret_cast<char*>(&energy),sizeof(double));
502 72 : pbc.setBox(box);
503 72 : }
504 :
505 53 : double Atoms::getKBoltzmann()const {
506 53 : if(naturalUnits || MDnaturalUnits) return 1.0;
507 46 : else return kBoltzmann/units.getEnergy();
508 : }
509 :
510 0 : double Atoms::getMDKBoltzmann()const {
511 0 : if(naturalUnits || MDnaturalUnits) return 1.0;
512 0 : else return kBoltzmann/MDUnits.getEnergy();
513 : }
514 :
515 0 : void Atoms::getLocalMasses(std::vector<double>& localMasses) {
516 0 : localMasses.resize(gatindex.size());
517 0 : for(unsigned i=0; i<gatindex.size(); i++) localMasses[i] = masses[gatindex[i]];
518 0 : }
519 :
520 0 : void Atoms::getLocalPositions(std::vector<Vector>& localPositions) {
521 0 : localPositions.resize(gatindex.size());
522 0 : mdatoms->getLocalPositions(localPositions);
523 0 : }
524 :
525 0 : void Atoms::getLocalForces(std::vector<Vector>& localForces) {
526 0 : localForces.resize(gatindex.size());
527 0 : for(unsigned i=0; i<gatindex.size(); i++) localForces[i] = forces[gatindex[i]];
528 0 : }
529 :
530 0 : void Atoms::getLocalMDForces(std::vector<Vector>& localForces) {
531 0 : localForces.resize(gatindex.size());
532 0 : for(unsigned i=0; i<gatindex.size(); i++) {
533 0 : localForces[i] = mdatoms->getMDforces(i);
534 : }
535 0 : }
536 :
537 2523 : }
|