ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/branches/new_design/OOPSE-3.0/src/brains/SimInfo.cpp
(Generate patch)

Comparing branches/new_design/OOPSE-3.0/src/brains/SimInfo.cpp (file contents):
Revision 1727 by tim, Thu Nov 11 16:41:58 2004 UTC vs.
Revision 1733 by tim, Fri Nov 12 06:19:04 2004 UTC

# Line 37 | Line 37 | SimInfo::SimInfo() : nAtoms_(0), nBonds_(0), nBends_(0
37  
38   namespace oopse {
39  
40 < SimInfo::SimInfo() : nAtoms_(0), nBonds_(0), nBends_(0), nTorsions_(0), nRigidBodies_(0),
41 <        nIntegrableObjects_(0), nCutoffGroups_(0), nConstraints_(0), sman_(NULL){
40 > SimInfo::SimInfo(const std::vector<std::pair<MoleculeStamp*, int> >& molStampPairs,
41 >                                ForceField* ff, Globals* globals) :
42 >                                forceField_(ff), globals_(globals), nAtoms_(0), nBonds_(0),
43 >                                nBends_(0), nTorsions_(0), nRigidBodies_(0), nIntegrableObjects_(0),
44 >                                nCutoffGroups_(0), nConstraints_(0), nZConstraint_(0), sman_(NULL) {
45  
46 +    std::vector<std::pair<MoleculeStamp*, int> >::iterator i;
47 +    int nCutoffAtoms; // number of atoms belong to cutoff groups
48 +    int ngroups;          //total cutoff groups defined in meta-data file
49 +    MoleculeStamp* molStamp;
50 +    int nMolWithSameStamp;
51 +    CutoffGroupStamp* cgStamp;
52 +    int nAtomsIngroups;
53 +    int nCutoffGroupsInStamp;    
54 +
55 +    nGlobalAtoms_ =  0;
56 +    ngroups = 0;
57 +    
58 +    for (i = molStampPairs.begin(); i !=molStampPairs.end(); ++i) {
59 +        molStamp = i->first;
60 +        nMolWithSameStamp = i->second;
61 +        
62 +        addMoleculeStamp(molStamp, nMolWithSameStamp);
63 +        
64 +        nGlobalAtoms_ += molStamp->getNAtoms() *nMolWithSameStamp;  
65 +        
66 +        nAtomsIngroups = 0;
67 +        nCutoffGroupsInStamp = molStamp->getNCutoffGroups();
68 +        
69 +        for (int j=0; j < nCutoffGroupsInStamp; j++) {
70 +            cgStamp = molStamp->getCutoffGroup(j);
71 +            nAtomsIngroups += cgStamp->getNMembers();
72 +        }
73 +
74 +        ngroups += *nMolWithSameStamp;
75 +        nCutoffAtoms += nAtomsIngroups * nMolWithSameStamp;                
76 +    }
77 +
78 +    //every free atom (atom does not belong to cutoff groups) is a cutoff group
79 +    //therefore the total number of cutoff groups in the system is equal to
80 +    //the total number of atoms minus number of atoms belong to cutoff group defined in meta-data
81 +    //file plus the number of cutoff groups defined in meta-data file
82 +    nGlobalCutoffGroups_ = nGlobalAtoms_ - nCutoffAtoms + ngroups;
83 +
84 +    //initialize globalGroupMembership_, every element of this array will be 0
85 +    globalGroupMembership_.insert(globalGroupMembership_.end(), nGlobalAtoms_, 0);
86 +
87 +    nGlobalMols_ = molStampIds_.size();
88 +
89 + #ifdef IS_MPI    
90 +    molToProcMap_.resize(nGlobalMols_);
91 + #endif
92 +    
93   }
94  
95   SimInfo::~SimInfo() {
96      //MemoryUtils::deleteVectorOfPointer(molecules_);
97 +
98 +    MemoryUtils::deleteVectorOfPointer(moleculeStamps_);
99 +    
100      delete sman_;
101 +    delete globals_;
102 +    delete forceField_;
103  
104   }
105  
106  
107   bool SimInfo::addMolecule(Molecule* mol) {
108      MoleculeIterator i;
109 <    i = std::find(molecules_.begin(), molecules_.end(), mol);
109 >
110 >    i = molecules_.find(mol->getGlobalIndex());
111      if (i != molecules_.end() ) {
112  
113          molecules_.insert(make_pair(mol->getGlobalIndex(), mol));
# Line 73 | Line 129 | bool SimInfo::removeMolecule(Molecule* mol) {
129  
130   bool SimInfo::removeMolecule(Molecule* mol) {
131      MoleculeIterator i;
132 <    i = std::find(molecules_.begin(), molecules_.end(), mol);
132 >    i = molecules_.find(mol->getGlobalIndex());
133  
134      if (i != molecules_.end() ) {
135  
136 +        assert(mol == i->second);
137 +        
138          nAtoms_ -= mol->getNAtoms();
139          nBonds_ -= mol->getNBonds();
140          nBends_ -= mol->getNBends();
# Line 145 | Line 203 | void SimInfo::calcNdf() {
203      ndf_ = ndf_local;
204   #endif
205  
206 <    // nZconstraints is global, as are the 3 COM translations for the
206 >    // nZconstraints_ is global, as are the 3 COM translations for the
207      // entire system:
208 <    ndf_ = ndf_ - 3 - nZconstraints;
208 >    ndf_ = ndf_ - 3 - nZconstraints_;
209  
210   }
211  
# Line 198 | Line 256 | void SimInfo::calcNdfTrans() {
256      ndfTrans_ = ndfTrans_local;
257   #endif
258  
259 <    ndfTrans_ = ndfTrans_ - 3 - nZconstraints;
259 >    ndfTrans_ = ndfTrans_ - 3 - nZconstraints_;
260  
261   }
262  
# Line 304 | Line 362 | void SimInfo::update() {
362  
363   void SimInfo::update() {
364  
365 + #ifdef IS_MPI
366 +    setupFortranParallel();
367 + #endif
368  
369 +    setupFortranSim();
370 +
371 +    calcNdf();
372 +    calcNdfRaw();
373 +    calcNdfTrans();
374 + }
375 +
376 + std::set<AtomType*> SimInfo::getUniqueAtomTypes() {
377 +    typename SimInfo::MoleculeIterator mi;
378 +    Molecule* mol;
379 +    typename Molecule::AtomIterator ai;
380 +    Atom* atom;
381 +    std::set<AtomType*> atomTypes;
382 +
383 +    for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {
384 +
385 +        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
386 +            atomTypes.insert(atom->getAtomType());
387 +        }
388 +        
389 +    }
390 +
391 +    return atomTypes;        
392 + }
393 +
394 + void SimInfo::setupSimType() {
395 +    std::set<AtomType*>::iterator i;
396 +    std::set<AtomType*> atomTypes;
397 +    atomTypes = getUniqueAtomTypes();
398      
399 +    int useLennardJones = 0;
400 +    int useElectrostatic = 0;
401 +    int useEAM = 0;
402 +    int useCharge = 0;
403 +    int useDirectional = 0;
404 +    int useDipole = 0;
405 +    int useGayBerne = 0;
406 +    int useSticky = 0;
407 +    int useShape = 0;
408 +    int useFLARB = 0; //it is not in AtomType yet
409 +    int useDirectionalAtom = 0;    
410 +    int useElectrostatics = 0;
411 +    //usePBC and useRF are from globals
412 +    bool usePBC = globals_->getPBC();
413 +    bool useRF = globals_->getUseRF();
414 +
415 +    //loop over all of the atom types
416 +    for (i = atomTypes.begin(); i != atomTypes.end(); ++i) {
417 +        useLennardJones |= i->isLennardJones();
418 +        useElectrostatic |= i->isElectrostatic();
419 +        useEAM |= i->isEAM();
420 +        useCharge |= i->isCharge();
421 +        useDirectional |= i->isDirectional();
422 +        useDipole |= i->isDipole();
423 +        useGayBerne |= i->isGayBerne();
424 +        useSticky |= i->isSticky();
425 +        useShape |= i->isShape();
426 +    }
427 +
428 +    if (useSticky || useDipole || useGayBerne || useShape) {
429 +        useDirectionalAtom = 1;
430 +    }
431 +
432 +    if (useCharge || useDipole) {
433 +        useElectrostatics = 1;
434 +    }
435 +
436 + #ifdef IS_MPI    
437 +    int temp;
438 +
439 +    temp = usePBC;
440 +    MPI_Allreduce(&temp, &usePBC, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
441 +
442 +    temp = useDirectionalAtom;
443 +    MPI_Allreduce(&temp, &useDirectionalAtom, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
444 +
445 +    temp = useLennardJones;
446 +    MPI_Allreduce(&temp, &useLennardJones, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
447 +
448 +    temp = useElectrostatics;
449 +    MPI_Allreduce(&temp, &useElectrostatics, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
450 +
451 +    temp = useCharge;
452 +    MPI_Allreduce(&temp, &useCharge, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
453 +
454 +    temp = useDipole;
455 +    MPI_Allreduce(&temp, &useDipole, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
456 +
457 +    temp = useSticky;
458 +    MPI_Allreduce(&temp, &useSticky, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
459 +
460 +    temp = useGayBerne;
461 +    MPI_Allreduce(&temp, &useGayBerne, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
462 +
463 +    temp = useEAM;
464 +    MPI_Allreduce(&temp, &useEAM, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
465 +
466 +    temp = useShape;
467 +    MPI_Allreduce(&temp, &useShape, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);  
468 +
469 +    temp = useFLARB;
470 +    MPI_Allreduce(&temp, &useFLARB, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
471 +
472 +    temp = useRF;
473 +    MPI_Allreduce(&temp, &useRF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
474 +    
475 + #endif
476 +
477 +    fInfo_.SIM_uses_PBC = usePBC;    
478 +    fInfo_.SIM_uses_DirectionalAtoms = useDirectionalAtom;
479 +    fInfo_.SIM_uses_LennardJones = useLennardJones;
480 +    fInfo_.SIM_uses_Electrostatics = useElectrostatics;    
481 +    fInfo_.SIM_uses_Charges = useCharge;
482 +    fInfo_.SIM_uses_Dipoles = useDipole;
483 +    fInfo_.SIM_uses_Sticky = useSticky;
484 +    fInfo_.SIM_uses_GayBerne = useGayBerne;
485 +    fInfo_.SIM_uses_EAM = useEAM;
486 +    fInfo_.SIM_uses_Shapes = useShape;
487 +    fInfo_.SIM_uses_FLARB = useFLARB;
488 +    fInfo_.SIM_uses_RF = useRF;
489 +
490 +    if( fInfo_.SIM_uses_Dipoles && fInfo_.SIM_uses_RF) {
491 +        fInfo_.dielect = dielectric;
492 +    } else {
493 +        fInfo_.dielect = 0.0;
494 +    }
495 +
496 + }
497 +
498 + void SimInfo::setupFortranSim() {
499 +    int isError;
500 +    int nExclude;
501 +    std::vector<int> fortranGlobalGroupMembership;
502 +    
503 +    nExclude = exclude_.getSize();
504 +    isError = 0;
505 +
506 +    //globalGroupMembership_ is filled by SimCreator    
507 +    for (int i = 0; i < nGlobalAtoms_; i++) {
508 +        fortranGlobalGroupMembership.push_back(globalGroupMembership_[i] + 1);
509 +    }
510 +
511 +    //calculate mass ratio of cutoff group
512 +    std::vector<double> mfact;
513 +    typename SimInfo::MoleculeIterator mi;
514 +    Molecule* mol;
515 +    typename Molecule::CutoffGroupIterator ci;
516 +    CutoffGroup* cg;
517 +    typename Molecule::AtomIterator ai;
518 +    Atom* atom;
519 +    double totalMass;
520 +
521 +    //to avoid memory reallocation, reserve enough space for mfact
522 +    mfact.reserve(getNCutoffGroups());
523 +    
524 +    for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {        
525 +        for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci)) {
526 +
527 +            totalMass = cg->getMass();
528 +            for(atom = cg->beginAtom(ai); atom != NULL; atom = cg->nextAtom(ai)) {
529 +                        mfact.push_back(atom->getMass()/totalMass);
530 +            }
531 +
532 +        }      
533 +    }
534 +
535 +    //fill ident array of local atoms (it is actually ident of AtomType, it is so confusing !!!)
536 +    std::vector<int> identArray;
537 +
538 +    //to avoid memory reallocation, reserve enough space identArray
539 +    identArray.reserve(getNAtoms());
540 +    
541 +    for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {        
542 +        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
543 +            identArray.push_back(atom->getIdent());
544 +        }
545 +    }    
546 +
547 +    //setup fortran simulation
548 +    //gloalExcludes and molMembershipArray should go away (They are never used)
549 +    //why the hell fortran need to know molecule?
550 +    //OOPSE = Object-Obfuscated Parallel Simulation Engine
551 +    
552 +    setFortranSim( &fInfo_, &nGlobalAtoms_, &nAtoms_, &identArray[0], &nExclude, exclude_->getExcludeList(),
553 +                  &nGlobalExcludes, globalExcludes, molMembershipArray,
554 +                  &mfact[0], &nCutoffGroups_, &fortranGlobalGroupMembership[0], &isError);
555 +
556 +    if( isError ){
557 +
558 +        sprintf( painCave.errMsg,
559 +                 "There was an error setting the simulation information in fortran.\n" );
560 +        painCave.isFatal = 1;
561 +        painCave.severity = OOPSE_ERROR;
562 +        simError();
563 +    }
564 +
565 + #ifdef IS_MPI
566 +    sprintf( checkPointMsg,
567 +       "succesfully sent the simulation information to fortran.\n");
568 +    MPIcheckPoint();
569 + #endif // is_mpi
570 + }
571 +
572 +
573 + #ifdef IS_MPI
574 + void SimInfo::setupFortranParallel() {
575 +    
576      //SimInfo is responsible for creating localToGlobalAtomIndex and localToGlobalGroupIndex
577      std::vector<int> localToGlobalAtomIndex(getNAtoms(), 0);
578      std::vector<int> localToGlobalCutoffGroupIndex;
# Line 332 | Line 599 | void SimInfo::update() {
599          
600      }
601  
602 <    //Setup Parallel Data and pass the index arrays to fortran
603 <    parallelData.nMolGlobal = getNMolGlobal();
604 <    parallelData.nMolLocal = ;
605 <    parallelData.nAtomsGlobal = ;
606 <    parallelData.nAtomsLocal = ;
607 <    parallelData.nGroupsGlobal = ;
608 <    parallelData.nGroupsLocal = ;
602 >    //fill up mpiSimData struct
603 >    parallelData.nMolGlobal = getNGlobalMolecules();
604 >    parallelData.nMolLocal = getNMolecules();
605 >    parallelData.nAtomsGlobal = getNGlobalAtoms();
606 >    parallelData.nAtomsLocal = getNAtoms();
607 >    parallelData.nGroupsGlobal = getNGlobalCutoffGroups();
608 >    parallelData.nGroupsLocal = getNCutoffGroups();
609      parallelData.myNode = worldRank;
610      MPI_Comm_size(MPI_COMM_WORLD, &(parallelData->nProcessors));
611 <    
612 <    setFsimParallel(parallelData,            &(parallelData->nAtomsLocal),
611 >
612 >    //pass mpiSimData struct and index arrays to fortran
613 >    setFsimParallel(parallelData, &(parallelData->nAtomsLocal),
614                      &localToGlobalAtomIndex[0],  &(parallelData->nGroupsLocal),
615                      &localToGlobalCutoffGroupIndex[0], &isError);
616  
# Line 359 | Line 627 | std::ostream& operator <<(ostream& o, SimInfo& info) {
627  
628   }
629  
630 + #endif
631 +
632 + void SimInfo::addProperty(GenericData* genData) {
633 +    properties_.addProperty(genData);  
634 + }
635 +
636 + void SimInfo::removeProperty(const std::string& propName) {
637 +    properties_.removeProperty(propName);  
638 + }
639 +
640 + void SimInfo::clearProperties() {
641 +    properties_.clearProperties();
642 + }
643 +
644 + std::vector<std::string> SimInfo::getPropertyNames() {
645 +    return properties_.getPropertyNames();  
646 + }
647 +      
648 + std::vector<GenericData*> SimInfo::getProperties() {
649 +    return properties_.getProperties();
650 + }
651 +
652 + GenericData* SimInfo::getPropertyByName(const std::string& propName) {
653 +    return properties_.getPropertyByName(propName);
654 + }
655 +
656 +
657   std::ostream& operator <<(ostream& o, SimInfo& info) {
658  
659      return o;

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines