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 1842 by tim, Fri Dec 3 20:30:07 2004 UTC

# Line 31 | Line 31
31   */
32  
33   #include <algorithm>
34 <
34 > #include <set>
35   #include "brains/SimInfo.hpp"
36 + #include "primitives/Molecule.hpp"
37 + #include "UseTheForce/doForces_interface.h"
38 + #include "UseTheForce/notifyCutoffs_interface.h"
39   #include "utils/MemoryUtils.hpp"
40 + #include "utils/simError.h"
41  
42   namespace oopse {
43  
44 < SimInfo::SimInfo() : nAtoms_(0), nBonds_(0), nBends_(0), nTorsions_(0), nRigidBodies_(0),
45 <        nIntegrableObjects_(0), nCutoffGroups_(0), nConstraints_(0), sman_(NULL){
44 > SimInfo::SimInfo(std::vector<std::pair<MoleculeStamp*, int> >& molStampPairs,
45 >                                ForceField* ff, Globals* simParams) :
46 >                                forceField_(ff), simParams_(simParams),
47 >                                ndf_(0), ndfRaw_(0), ndfTrans_(0), nZconstraint_(0),
48 >                                nGlobalMols_(0), nGlobalAtoms_(0), nGlobalCutoffGroups_(0),
49 >                                nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0),
50 >                                nAtoms_(0), nBonds_(0),  nBends_(0), nTorsions_(0), nRigidBodies_(0),
51 >                                nIntegrableObjects_(0),  nCutoffGroups_(0), nConstraints_(0),
52 >                                sman_(NULL), fortranInitialized_(false) {
53 >
54 >            
55 >    std::vector<std::pair<MoleculeStamp*, int> >::iterator i;
56 >    MoleculeStamp* molStamp;
57 >    int nMolWithSameStamp;
58 >    int nCutoffAtoms = 0; // number of atoms belong to cutoff groups
59 >    int nGroups = 0;          //total cutoff groups defined in meta-data file
60 >    CutoffGroupStamp* cgStamp;    
61 >    RigidBodyStamp* rbStamp;
62 >    int nRigidAtoms = 0;
63 >    
64 >    for (i = molStampPairs.begin(); i !=molStampPairs.end(); ++i) {
65 >        molStamp = i->first;
66 >        nMolWithSameStamp = i->second;
67 >        
68 >        addMoleculeStamp(molStamp, nMolWithSameStamp);
69 >
70 >        //calculate atoms in molecules
71 >        nGlobalAtoms_ += molStamp->getNAtoms() *nMolWithSameStamp;  
72 >
73 >
74 >        //calculate atoms in cutoff groups
75 >        int nAtomsInGroups = 0;
76 >        int nCutoffGroupsInStamp = molStamp->getNCutoffGroups();
77 >        
78 >        for (int j=0; j < nCutoffGroupsInStamp; j++) {
79 >            cgStamp = molStamp->getCutoffGroup(j);
80 >            nAtomsInGroups += cgStamp->getNMembers();
81 >        }
82 >
83 >        nGroups += nCutoffGroupsInStamp * nMolWithSameStamp;
84 >        nCutoffAtoms += nAtomsInGroups * nMolWithSameStamp;            
85 >
86 >        //calculate atoms in rigid bodies
87 >        int nAtomsInRigidBodies = 0;
88 >        int nRigidBodiesInStamp = molStamp->getNCutoffGroups();
89 >        
90 >        for (int j=0; j < nRigidBodiesInStamp; j++) {
91 >            rbStamp = molStamp->getRigidBody(j);
92 >            nAtomsInRigidBodies += rbStamp->getNMembers();
93 >        }
94 >
95 >        nGlobalRigidBodies_ += nRigidBodiesInStamp * nMolWithSameStamp;
96 >        nRigidAtoms += nAtomsInRigidBodies * nMolWithSameStamp;            
97 >        
98 >    }
99 >
100 >    //every free atom (atom does not belong to cutoff groups) is a cutoff group
101 >    //therefore the total number of cutoff groups in the system is equal to
102 >    //the total number of atoms minus number of atoms belong to cutoff group defined in meta-data
103 >    //file plus the number of cutoff groups defined in meta-data file
104 >    nGlobalCutoffGroups_ = nGlobalAtoms_ - nCutoffAtoms + nGroups;
105  
106 +    //every free atom (atom does not belong to rigid bodies) is an integrable object
107 +    //therefore the total number of  integrable objects in the system is equal to
108 +    //the total number of atoms minus number of atoms belong to  rigid body defined in meta-data
109 +    //file plus the number of  rigid bodies defined in meta-data file
110 +    nGlobalIntegrableObjects_ = nGlobalAtoms_ - nRigidAtoms + nGlobalRigidBodies_;
111 +
112 +    nGlobalMols_ = molStampIds_.size();
113 +
114 + #ifdef IS_MPI    
115 +    molToProcMap_.resize(nGlobalMols_);
116 + #endif
117 +    
118   }
119  
120   SimInfo::~SimInfo() {
121      //MemoryUtils::deleteVectorOfPointer(molecules_);
122 +
123 +    MemoryUtils::deleteVectorOfPointer(moleculeStamps_);
124 +    
125      delete sman_;
126 +    delete simParams_;
127 +    delete forceField_;
128  
129   }
130  
131  
132   bool SimInfo::addMolecule(Molecule* mol) {
133      MoleculeIterator i;
54    i = std::find(molecules_.begin(), molecules_.end(), mol);
55    if (i != molecules_.end() ) {
134  
135 <        molecules_.insert(make_pair(mol->getGlobalIndex(), mol));
135 >    i = molecules_.find(mol->getGlobalIndex());
136 >    if (i == molecules_.end() ) {
137 >
138 >        molecules_.insert(std::make_pair(mol->getGlobalIndex(), mol));
139          
140          nAtoms_ += mol->getNAtoms();
141          nBonds_ += mol->getNBonds();
# Line 73 | Line 154 | bool SimInfo::removeMolecule(Molecule* mol) {
154  
155   bool SimInfo::removeMolecule(Molecule* mol) {
156      MoleculeIterator i;
157 <    i = std::find(molecules_.begin(), molecules_.end(), mol);
157 >    i = molecules_.find(mol->getGlobalIndex());
158  
159      if (i != molecules_.end() ) {
160  
161 +        assert(mol == i->second);
162 +        
163          nAtoms_ -= mol->getNAtoms();
164          nBonds_ -= mol->getNBonds();
165          nBends_ -= mol->getNBends();
# Line 101 | Line 184 | Molecule* SimInfo::beginMolecule(MoleculeIterator& i)
184          
185   Molecule* SimInfo::beginMolecule(MoleculeIterator& i) {
186      i = molecules_.begin();
187 <    return i == molecules_.end() ? NULL : *i;
187 >    return i == molecules_.end() ? NULL : i->second;
188   }    
189  
190   Molecule* SimInfo::nextMolecule(MoleculeIterator& i) {
191      ++i;
192 <    return i == molecules_.end() ? NULL : *i;    
192 >    return i == molecules_.end() ? NULL : i->second;    
193   }
194  
195  
# Line 145 | Line 228 | void SimInfo::calcNdf() {
228      ndf_ = ndf_local;
229   #endif
230  
231 <    // nZconstraints is global, as are the 3 COM translations for the
231 >    // nZconstraints_ is global, as are the 3 COM translations for the
232      // entire system:
233 <    ndf_ = ndf_ - 3 - nZconstraints;
233 >    ndf_ = ndf_ - 3 - nZconstraint_;
234  
235   }
236  
# Line 198 | Line 281 | void SimInfo::calcNdfTrans() {
281      ndfTrans_ = ndfTrans_local;
282   #endif
283  
284 <    ndfTrans_ = ndfTrans_ - 3 - nZconstraints;
284 >    ndfTrans_ = ndfTrans_ - 3 - nZconstraint_;
285  
286   }
287  
# Line 230 | Line 313 | void SimInfo::addExcludePairs(Molecule* mol) {
313          exclude_.addPair(b, c);        
314      }
315  
316 <    for (torsion= mol->beginTorsion(torsionIter); torsion != NULL; torsion = mol->nextBond(torsionIter)) {
316 >    for (torsion= mol->beginTorsion(torsionIter); torsion != NULL; torsion = mol->nextTorsion(torsionIter)) {
317          a = torsion->getAtomA()->getGlobalIndex();
318          b = torsion->getAtomB()->getGlobalIndex();        
319          c = torsion->getAtomC()->getGlobalIndex();        
# Line 275 | Line 358 | void SimInfo::removeExcludePairs(Molecule* mol) {
358          exclude_.removePair(b, c);        
359      }
360  
361 <    for (torsion= mol->beginTorsion(torsionIter); torsion != NULL; torsion = mol->nextBond(torsionIter)) {
361 >    for (torsion= mol->beginTorsion(torsionIter); torsion != NULL; torsion = mol->nextTorsion(torsionIter)) {
362          a = torsion->getAtomA()->getGlobalIndex();
363          b = torsion->getAtomB()->getGlobalIndex();        
364          c = torsion->getAtomC()->getGlobalIndex();        
# Line 299 | Line 382 | void SimInfo::addMoleculeStamp(MoleculeStamp* molStamp
382      curStampId = molStampIds_.size();
383  
384      moleculeStamps_.push_back(molStamp);
385 <    molStampIds_.insert(molStampIds_.end(), nmol, curStampId)
385 >    molStampIds_.insert(molStampIds_.end(), nmol, curStampId);
386   }
387  
388   void SimInfo::update() {
389  
390 +    setupSimType();
391  
392 + #ifdef IS_MPI
393 +    setupFortranParallel();
394 + #endif
395 +
396 +    setupFortranSim();
397 +
398 +    //setup fortran force field
399 +    /** @deprecate */    
400 +    int isError = 0;
401 +    initFortranFF( &fInfo_.SIM_uses_RF , &isError );
402 +    if(isError){
403 +        sprintf( painCave.errMsg,
404 +         "ForceField error: There was an error initializing the forceField in fortran.\n" );
405 +        painCave.isFatal = 1;
406 +        simError();
407 +    }
408 +  
409      
410 +    setupCutoff();
411 +
412 +    calcNdf();
413 +    calcNdfRaw();
414 +    calcNdfTrans();
415 +
416 +    fortranInitialized_ = true;
417 + }
418 +
419 + std::set<AtomType*> SimInfo::getUniqueAtomTypes() {
420 +    SimInfo::MoleculeIterator mi;
421 +    Molecule* mol;
422 +    Molecule::AtomIterator ai;
423 +    Atom* atom;
424 +    std::set<AtomType*> atomTypes;
425 +
426 +    for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {
427 +
428 +        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
429 +            atomTypes.insert(atom->getAtomType());
430 +        }
431 +        
432 +    }
433 +
434 +    return atomTypes;        
435 + }
436 +
437 + void SimInfo::setupSimType() {
438 +    std::set<AtomType*>::iterator i;
439 +    std::set<AtomType*> atomTypes;
440 +    atomTypes = getUniqueAtomTypes();
441 +    
442 +    int useLennardJones = 0;
443 +    int useElectrostatic = 0;
444 +    int useEAM = 0;
445 +    int useCharge = 0;
446 +    int useDirectional = 0;
447 +    int useDipole = 0;
448 +    int useGayBerne = 0;
449 +    int useSticky = 0;
450 +    int useShape = 0;
451 +    int useFLARB = 0; //it is not in AtomType yet
452 +    int useDirectionalAtom = 0;    
453 +    int useElectrostatics = 0;
454 +    //usePBC and useRF are from simParams
455 +    bool usePBC = simParams_->getPBC();
456 +    bool useRF = simParams_->getUseRF();
457 +
458 +    //loop over all of the atom types
459 +    for (i = atomTypes.begin(); i != atomTypes.end(); ++i) {
460 +        useLennardJones |= (*i)->isLennardJones();
461 +        useElectrostatic |= (*i)->isElectrostatic();
462 +        useEAM |= (*i)->isEAM();
463 +        useCharge |= (*i)->isCharge();
464 +        useDirectional |= (*i)->isDirectional();
465 +        useDipole |= (*i)->isDipole();
466 +        useGayBerne |= (*i)->isGayBerne();
467 +        useSticky |= (*i)->isSticky();
468 +        useShape |= (*i)->isShape();
469 +    }
470 +
471 +    if (useSticky || useDipole || useGayBerne || useShape) {
472 +        useDirectionalAtom = 1;
473 +    }
474 +
475 +    if (useCharge || useDipole) {
476 +        useElectrostatics = 1;
477 +    }
478 +
479 + #ifdef IS_MPI    
480 +    int temp;
481 +
482 +    temp = usePBC;
483 +    MPI_Allreduce(&temp, &usePBC, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
484 +
485 +    temp = useDirectionalAtom;
486 +    MPI_Allreduce(&temp, &useDirectionalAtom, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
487 +
488 +    temp = useLennardJones;
489 +    MPI_Allreduce(&temp, &useLennardJones, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
490 +
491 +    temp = useElectrostatics;
492 +    MPI_Allreduce(&temp, &useElectrostatics, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
493 +
494 +    temp = useCharge;
495 +    MPI_Allreduce(&temp, &useCharge, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
496 +
497 +    temp = useDipole;
498 +    MPI_Allreduce(&temp, &useDipole, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
499 +
500 +    temp = useSticky;
501 +    MPI_Allreduce(&temp, &useSticky, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
502 +
503 +    temp = useGayBerne;
504 +    MPI_Allreduce(&temp, &useGayBerne, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
505 +
506 +    temp = useEAM;
507 +    MPI_Allreduce(&temp, &useEAM, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
508 +
509 +    temp = useShape;
510 +    MPI_Allreduce(&temp, &useShape, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);  
511 +
512 +    temp = useFLARB;
513 +    MPI_Allreduce(&temp, &useFLARB, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
514 +
515 +    temp = useRF;
516 +    MPI_Allreduce(&temp, &useRF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
517 +    
518 + #endif
519 +
520 +    fInfo_.SIM_uses_PBC = usePBC;    
521 +    fInfo_.SIM_uses_DirectionalAtoms = useDirectionalAtom;
522 +    fInfo_.SIM_uses_LennardJones = useLennardJones;
523 +    fInfo_.SIM_uses_Electrostatics = useElectrostatics;    
524 +    fInfo_.SIM_uses_Charges = useCharge;
525 +    fInfo_.SIM_uses_Dipoles = useDipole;
526 +    fInfo_.SIM_uses_Sticky = useSticky;
527 +    fInfo_.SIM_uses_GayBerne = useGayBerne;
528 +    fInfo_.SIM_uses_EAM = useEAM;
529 +    fInfo_.SIM_uses_Shapes = useShape;
530 +    fInfo_.SIM_uses_FLARB = useFLARB;
531 +    fInfo_.SIM_uses_RF = useRF;
532 +
533 +    if( fInfo_.SIM_uses_Dipoles && fInfo_.SIM_uses_RF) {
534 +
535 +        if (simParams_->haveDielectric()) {
536 +            fInfo_.dielect = simParams_->getDielectric();
537 +        } else {
538 +            sprintf(painCave.errMsg,
539 +                    "SimSetup Error: No Dielectric constant was set.\n"
540 +                    "\tYou are trying to use Reaction Field without"
541 +                    "\tsetting a dielectric constant!\n");
542 +            painCave.isFatal = 1;
543 +            simError();
544 +        }
545 +        
546 +    } else {
547 +        fInfo_.dielect = 0.0;
548 +    }
549 +
550 + }
551 +
552 + void SimInfo::setupFortranSim() {
553 +    int isError;
554 +    int nExclude;
555 +    std::vector<int> fortranGlobalGroupMembership;
556 +    
557 +    nExclude = exclude_.getSize();
558 +    isError = 0;
559 +
560 +    //globalGroupMembership_ is filled by SimCreator    
561 +    for (int i = 0; i < nGlobalAtoms_; i++) {
562 +        fortranGlobalGroupMembership.push_back(globalGroupMembership_[i] + 1);
563 +    }
564 +
565 +    //calculate mass ratio of cutoff group
566 +    std::vector<double> mfact;
567 +    SimInfo::MoleculeIterator mi;
568 +    Molecule* mol;
569 +    Molecule::CutoffGroupIterator ci;
570 +    CutoffGroup* cg;
571 +    Molecule::AtomIterator ai;
572 +    Atom* atom;
573 +    double totalMass;
574 +
575 +    //to avoid memory reallocation, reserve enough space for mfact
576 +    mfact.reserve(getNCutoffGroups());
577 +    
578 +    for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {        
579 +        for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci)) {
580 +
581 +            totalMass = cg->getMass();
582 +            for(atom = cg->beginAtom(ai); atom != NULL; atom = cg->nextAtom(ai)) {
583 +                        mfact.push_back(atom->getMass()/totalMass);
584 +            }
585 +
586 +        }      
587 +    }
588 +
589 +    //fill ident array of local atoms (it is actually ident of AtomType, it is so confusing !!!)
590 +    std::vector<int> identArray;
591 +
592 +    //to avoid memory reallocation, reserve enough space identArray
593 +    identArray.reserve(getNAtoms());
594 +    
595 +    for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {        
596 +        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
597 +            identArray.push_back(atom->getIdent());
598 +        }
599 +    }    
600 +
601 +    //fill molMembershipArray
602 +    //molMembershipArray is filled by SimCreator    
603 +    std::vector<int> molMembershipArray(nGlobalAtoms_);
604 +    for (int i = 0; i < nGlobalAtoms_; i++) {
605 +        molMembershipArray.push_back(globalMolMembership_[i] + 1);
606 +    }
607 +    
608 +    //setup fortran simulation
609 +    //gloalExcludes and molMembershipArray should go away (They are never used)
610 +    //why the hell fortran need to know molecule?
611 +    //OOPSE = Object-Obfuscated Parallel Simulation Engine
612 +    int nGlobalExcludes = 0;
613 +    int* globalExcludes = NULL;
614 +    int* excludeList = exclude_.getExcludeList();
615 +    setFortranSim( &fInfo_, &nGlobalAtoms_, &nAtoms_, &identArray[0], &nExclude, excludeList ,
616 +                  &nGlobalExcludes, globalExcludes, &molMembershipArray[0],
617 +                  &mfact[0], &nCutoffGroups_, &fortranGlobalGroupMembership[0], &isError);
618 +
619 +    if( isError ){
620 +
621 +        sprintf( painCave.errMsg,
622 +                 "There was an error setting the simulation information in fortran.\n" );
623 +        painCave.isFatal = 1;
624 +        painCave.severity = OOPSE_ERROR;
625 +        simError();
626 +    }
627 +
628 + #ifdef IS_MPI
629 +    sprintf( checkPointMsg,
630 +       "succesfully sent the simulation information to fortran.\n");
631 +    MPIcheckPoint();
632 + #endif // is_mpi
633 + }
634 +
635 +
636 + #ifdef IS_MPI
637 + void SimInfo::setupFortranParallel() {
638 +    
639      //SimInfo is responsible for creating localToGlobalAtomIndex and localToGlobalGroupIndex
640      std::vector<int> localToGlobalAtomIndex(getNAtoms(), 0);
641      std::vector<int> localToGlobalCutoffGroupIndex;
642 <    typename SimInfo::MoleculeIterator mi;
643 <    typename Molecule::AtomIterator ai;
644 <    typename Molecule::CutoffGroupIterator ci;
642 >    SimInfo::MoleculeIterator mi;
643 >    Molecule::AtomIterator ai;
644 >    Molecule::CutoffGroupIterator ci;
645      Molecule* mol;
646      Atom* atom;
647      CutoffGroup* cg;
# Line 332 | Line 662 | void SimInfo::update() {
662          
663      }
664  
665 <    //Setup Parallel Data and pass the index arrays to fortran
666 <    parallelData.nMolGlobal = getNMolGlobal();
667 <    parallelData.nMolLocal = ;
668 <    parallelData.nAtomsGlobal = ;
669 <    parallelData.nAtomsLocal = ;
670 <    parallelData.nGroupsGlobal = ;
671 <    parallelData.nGroupsLocal = ;
665 >    //fill up mpiSimData struct
666 >    parallelData.nMolGlobal = getNGlobalMolecules();
667 >    parallelData.nMolLocal = getNMolecules();
668 >    parallelData.nAtomsGlobal = getNGlobalAtoms();
669 >    parallelData.nAtomsLocal = getNAtoms();
670 >    parallelData.nGroupsGlobal = getNGlobalCutoffGroups();
671 >    parallelData.nGroupsLocal = getNCutoffGroups();
672      parallelData.myNode = worldRank;
673      MPI_Comm_size(MPI_COMM_WORLD, &(parallelData->nProcessors));
674 <    
675 <    setFsimParallel(parallelData,            &(parallelData->nAtomsLocal),
674 >
675 >    //pass mpiSimData struct and index arrays to fortran
676 >    setFsimParallel(parallelData, &(parallelData->nAtomsLocal),
677                      &localToGlobalAtomIndex[0],  &(parallelData->nGroupsLocal),
678                      &localToGlobalCutoffGroupIndex[0], &isError);
679  
# Line 359 | Line 690 | std::ostream& operator <<(ostream& o, SimInfo& info) {
690  
691   }
692  
693 < std::ostream& operator <<(ostream& o, SimInfo& info) {
693 > #endif
694  
695 + double SimInfo::calcMaxCutoffRadius() {
696 +
697 +
698 +    std::set<AtomType*> atomTypes;
699 +    std::set<AtomType*>::iterator i;
700 +    std::vector<double> cutoffRadius;
701 +
702 +    //get the unique atom types
703 +    atomTypes = getUniqueAtomTypes();
704 +
705 +    //query the max cutoff radius among these atom types
706 +    for (i = atomTypes.begin(); i != atomTypes.end(); ++i) {
707 +        cutoffRadius.push_back(forceField_->getRcutFromAtomType(*i));
708 +    }
709 +
710 +    double maxCutoffRadius = *(std::max_element(cutoffRadius.begin(), cutoffRadius.end()));
711 + #ifdef IS_MPI
712 +    //pick the max cutoff radius among the processors
713 + #endif
714 +
715 +    return maxCutoffRadius;
716 + }
717 +
718 + void SimInfo::setupCutoff() {
719 +    double rcut_;  //cutoff radius
720 +    double rsw_; //switching radius
721 +    
722 +    if (fInfo_.SIM_uses_Charges | fInfo_.SIM_uses_Dipoles | fInfo_.SIM_uses_RF) {
723 +        
724 +        if (!simParams_->haveRcut()){
725 +            sprintf(painCave.errMsg,
726 +                "SimCreator Warning: No value was set for the cutoffRadius.\n"
727 +                "\tOOPSE will use a default value of 15.0 angstroms"
728 +                "\tfor the cutoffRadius.\n");
729 +            painCave.isFatal = 0;
730 +            simError();
731 +            rcut_ = 15.0;
732 +        } else{
733 +            rcut_ = simParams_->getRcut();
734 +        }
735 +
736 +        if (!simParams_->haveRsw()){
737 +            sprintf(painCave.errMsg,
738 +                "SimCreator Warning: No value was set for switchingRadius.\n"
739 +                "\tOOPSE will use a default value of\n"
740 +                "\t0.95 * cutoffRadius for the switchingRadius\n");
741 +            painCave.isFatal = 0;
742 +            simError();
743 +            rsw_ = 0.95 * rcut_;
744 +        } else{
745 +            rsw_ = simParams_->getRsw();
746 +        }
747 +
748 +    } else {
749 +        // if charge, dipole or reaction field is not used and the cutofff radius is not specified in
750 +        //meta-data file, the maximum cutoff radius calculated from forcefiled will be used
751 +        
752 +        if (simParams_->haveRcut()) {
753 +            rcut_ = simParams_->getRcut();
754 +        } else {
755 +            //set cutoff radius to the maximum cutoff radius based on atom types in the whole system
756 +            rcut_ = calcMaxCutoffRadius();
757 +        }
758 +
759 +        if (simParams_->haveRsw()) {
760 +            rsw_  = simParams_->getRsw();
761 +        } else {
762 +            rsw_ = rcut_;
763 +        }
764 +    
765 +    }
766 +        
767 +    double rnblist = rcut_ + 1; // skin of neighbor list
768 +
769 +    //Pass these cutoff radius etc. to fortran. This function should be called once and only once
770 +    notifyFortranCutoffs(&rcut_, &rsw_, &rnblist);
771 + }
772 +
773 + void SimInfo::addProperty(GenericData* genData) {
774 +    properties_.addProperty(genData);  
775 + }
776 +
777 + void SimInfo::removeProperty(const std::string& propName) {
778 +    properties_.removeProperty(propName);  
779 + }
780 +
781 + void SimInfo::clearProperties() {
782 +    properties_.clearProperties();
783 + }
784 +
785 + std::vector<std::string> SimInfo::getPropertyNames() {
786 +    return properties_.getPropertyNames();  
787 + }
788 +      
789 + std::vector<GenericData*> SimInfo::getProperties() {
790 +    return properties_.getProperties();
791 + }
792 +
793 + GenericData* SimInfo::getPropertyByName(const std::string& propName) {
794 +    return properties_.getPropertyByName(propName);
795 + }
796 +
797 + void SimInfo::setSnapshotManager(SnapshotManager* sman) {
798 +    sman_ = sman;
799 +
800 +    Molecule* mol;
801 +    RigidBody* rb;
802 +    Atom* atom;
803 +    SimInfo::MoleculeIterator mi;
804 +    Molecule::RigidBodyIterator rbIter;
805 +    Molecule::AtomIterator atomIter;;
806 +
807 +    for (mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {
808 +        
809 +        for (atom = mol->beginAtom(atomIter); atom != NULL; atom = mol->nextAtom(atomIter)) {
810 +            atom->setSnapshotManager(sman_);
811 +        }
812 +        
813 +        for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) {
814 +            rb->setSnapshotManager(sman_);
815 +        }
816 +    }    
817 +    
818 + }
819 +
820 + std::ostream& operator <<(std::ostream& o, SimInfo& info) {
821 +
822      return o;
823   }
824  
825   }//end namespace oopse
826 +

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines