# | Line 42 | Line 42 | SimInfo::SimInfo(){ | |
---|---|---|
42 | thermalTime = 0.0; | |
43 | currentTime = 0.0; | |
44 | rCut = 0.0; | |
45 | < | ecr = 0.0; |
46 | < | est = 0.0; |
45 | > | rSw = 0.0; |
46 | ||
47 | haveRcut = 0; | |
48 | < | haveEcr = 0; |
48 | > | haveRsw = 0; |
49 | boxIsInit = 0; | |
50 | ||
51 | resetTime = 1e99; | |
# | Line 63 | Line 62 | SimInfo::SimInfo(){ | |
62 | useReactionField = 0; | |
63 | useGB = 0; | |
64 | useEAM = 0; | |
65 | < | useMolecularCutoffs = 0; |
65 | > | useSolidThermInt = 0; |
66 | > | useLiquidThermInt = 0; |
67 | ||
68 | + | haveCutoffGroups = false; |
69 | + | |
70 | excludes = Exclude::Instance(); | |
71 | ||
72 | myConfiguration = new SimState(); | |
# | Line 189 | Line 191 | void SimInfo::calcHmatInv( void ) { | |
191 | ||
192 | if( oldOrtho != orthoRhombic ){ | |
193 | ||
194 | < | if( orthoRhombic ){ |
194 | > | if( orthoRhombic ) { |
195 | sprintf( painCave.errMsg, | |
196 | < | "OOPSE is switching from the default Non-Orthorhombic\n" |
196 | > | "\n\tOOPSE is switching from the default Non-Orthorhombic\n" |
197 | "\tto the faster Orthorhombic periodic boundary computations.\n" | |
198 | "\tThis is usually a good thing, but if you wan't the\n" | |
199 | "\tNon-Orthorhombic computations, make the orthoBoxTolerance\n" | |
200 | "\tvariable ( currently set to %G ) smaller.\n", | |
201 | orthoTolerance); | |
202 | + | painCave.severity = OOPSE_INFO; |
203 | simError(); | |
204 | } | |
205 | else { | |
206 | sprintf( painCave.errMsg, | |
207 | < | "OOPSE is switching from the faster Orthorhombic to the more\n" |
207 | > | "\n\tOOPSE is switching from the faster Orthorhombic to the more\n" |
208 | "\tflexible Non-Orthorhombic periodic boundary computations.\n" | |
209 | "\tThis is usually because the box has deformed under\n" | |
210 | "\tNPTf integration. If you wan't to live on the edge with\n" | |
211 | "\tthe Orthorhombic computations, make the orthoBoxTolerance\n" | |
212 | "\tvariable ( currently set to %G ) larger.\n", | |
213 | orthoTolerance); | |
214 | + | painCave.severity = OOPSE_WARNING; |
215 | simError(); | |
216 | } | |
217 | } | |
# | Line 444 | Line 448 | void SimInfo::refreshSim(){ | |
448 | excl = excludes->getFortranArray(); | |
449 | ||
450 | #ifdef IS_MPI | |
451 | < | n_global = mpiSim->getTotAtoms(); |
451 | > | n_global = mpiSim->getNAtomsGlobal(); |
452 | #else | |
453 | n_global = n_atoms; | |
454 | #endif | |
455 | ||
456 | isError = 0; | |
457 | ||
458 | < | getFortranGroupArray(this, mfact, ngroup, groupList, groupStart); |
459 | < | |
458 | > | getFortranGroupArrays(this, FglobalGroupMembership, mfact); |
459 | > | //it may not be a good idea to pass the address of first element in vector |
460 | > | //since c++ standard does not require vector to be stored continuously in meomory |
461 | > | //Most of the compilers will organize the memory of vector continuously |
462 | setFsimulation( &fInfo, &n_global, &n_atoms, identArray, &n_exclude, excl, | |
463 | &nGlobalExcludes, globalExcludes, molMembershipArray, | |
464 | < | &mfact[0], &ngroup, &groupList[0], &groupStart[0], &isError); |
465 | < | |
464 | > | &mfact[0], &ngroup, &FglobalGroupMembership[0], &isError); |
465 | > | |
466 | if( isError ){ | |
467 | ||
468 | sprintf( painCave.errMsg, | |
# | Line 480 | Line 486 | void SimInfo::setDefaultRcut( double theRcut ){ | |
486 | ||
487 | haveRcut = 1; | |
488 | rCut = theRcut; | |
489 | + | rList = rCut + 1.0; |
490 | ||
491 | < | ( rCut > ecr )? rList = rCut + 1.0: rList = ecr + 1.0; |
485 | < | |
486 | < | notifyFortranCutOffs( &rCut, &rList, &ecr, &est ); |
491 | > | notifyFortranCutOffs( &rCut, &rSw, &rList ); |
492 | } | |
493 | ||
494 | < | void SimInfo::setDefaultEcr( double theEcr ){ |
490 | < | |
491 | < | haveEcr = 1; |
492 | < | ecr = theEcr; |
493 | < | |
494 | < | ( rCut > ecr )? rList = rCut + 1.0: rList = ecr + 1.0; |
495 | < | |
496 | < | notifyFortranCutOffs( &rCut, &rList, &ecr, &est ); |
497 | < | } |
498 | < | |
499 | < | void SimInfo::setDefaultEcr( double theEcr, double theEst ){ |
494 | > | void SimInfo::setDefaultRcut( double theRcut, double theRsw ){ |
495 | ||
496 | < | est = theEst; |
497 | < | setDefaultEcr( theEcr ); |
496 | > | rSw = theRsw; |
497 | > | setDefaultRcut( theRcut ); |
498 | } | |
499 | ||
500 | ||
# | Line 511 | Line 506 | void SimInfo::checkCutOffs( void ){ | |
506 | ||
507 | if( rCut > maxCutoff ){ | |
508 | sprintf( painCave.errMsg, | |
509 | < | "LJrcut is too large for the current periodic box.\n" |
510 | < | "\tCurrent Value of LJrcut = %G at time %G\n " |
509 | > | "\n\tcutoffRadius is too large for the current periodic box.\n" |
510 | > | "\tCurrent Value of cutoffRadius = %G at time %G\n " |
511 | "\tThis is larger than half of at least one of the\n" | |
512 | "\tperiodic box vectors. Right now, the Box matrix is:\n" | |
513 | "\n" | |
# | Line 523 | Line 518 | void SimInfo::checkCutOffs( void ){ | |
518 | Hmat[0][0], Hmat[0][1], Hmat[0][2], | |
519 | Hmat[1][0], Hmat[1][1], Hmat[1][2], | |
520 | Hmat[2][0], Hmat[2][1], Hmat[2][2]); | |
521 | + | painCave.severity = OOPSE_ERROR; |
522 | painCave.isFatal = 1; | |
523 | simError(); | |
524 | < | } |
529 | < | |
530 | < | if( haveEcr ){ |
531 | < | if( ecr > maxCutoff ){ |
532 | < | sprintf( painCave.errMsg, |
533 | < | "electrostaticCutoffRadius is too large for the current\n" |
534 | < | "\tperiodic box.\n\n" |
535 | < | "\tCurrent Value of ECR = %G at time %G\n " |
536 | < | "\tThis is larger than half of at least one of the\n" |
537 | < | "\tperiodic box vectors. Right now, the Box matrix is:\n" |
538 | < | "\n" |
539 | < | "\t[ %G %G %G ]\n" |
540 | < | "\t[ %G %G %G ]\n" |
541 | < | "\t[ %G %G %G ]\n", |
542 | < | ecr, currentTime, |
543 | < | Hmat[0][0], Hmat[0][1], Hmat[0][2], |
544 | < | Hmat[1][0], Hmat[1][1], Hmat[1][2], |
545 | < | Hmat[2][0], Hmat[2][1], Hmat[2][2]); |
546 | < | painCave.isFatal = 1; |
547 | < | simError(); |
548 | < | } |
549 | < | } |
524 | > | } |
525 | } else { | |
526 | // initialize this stuff before using it, OK? | |
527 | sprintf( painCave.errMsg, | |
528 | < | "Trying to check cutoffs without a box.\n" |
528 | > | "\n\tTrying to check cutoffs without a box.\n" |
529 | "\tOOPSE should have better programmers than that.\n" ); | |
530 | + | painCave.severity = OOPSE_ERROR; |
531 | painCave.isFatal = 1; | |
532 | simError(); | |
533 | } | |
# | Line 595 | Line 571 | GenericData* SimInfo::getProperty(const string& propNa | |
571 | } | |
572 | ||
573 | ||
574 | < | void getFortranGroupArray(SimInfo* info, vector<double>& mfact, int& ngroup, |
575 | < | vector<int>& groupList, vector<int>& groupStart){ |
576 | < | Molecule* mol; |
574 | > | void SimInfo::getFortranGroupArrays(SimInfo* info, |
575 | > | vector<int>& FglobalGroupMembership, |
576 | > | vector<double>& mfact){ |
577 | > | |
578 | > | Molecule* myMols; |
579 | Atom** myAtoms; | |
580 | int numAtom; | |
603 | – | int curIndex; |
581 | double mtot; | |
582 | < | |
582 | > | int numMol; |
583 | > | int numCutoffGroups; |
584 | > | CutoffGroup* myCutoffGroup; |
585 | > | vector<CutoffGroup*>::iterator iterCutoff; |
586 | > | Atom* cutoffAtom; |
587 | > | vector<Atom*>::iterator iterAtom; |
588 | > | int atomIndex; |
589 | > | double totalMass; |
590 | > | |
591 | mfact.clear(); | |
592 | < | groupList.clear(); |
608 | < | groupStart.clear(); |
592 | > | FglobalGroupMembership.clear(); |
593 | ||
594 | < | //Be careful, fortran array begin at 1 |
595 | < | curIndex = 1; |
612 | < | |
613 | < | if(info->useMolecularCutoffs){ |
614 | < | |
594 | > | |
595 | > | // Fix the silly fortran indexing problem |
596 | #ifdef IS_MPI | |
597 | < | ngroup = mpiSim->getMyNMol(); |
597 | > | numAtom = mpiSim->getNAtomsGlobal(); |
598 | #else | |
599 | < | ngroup = info->n_mol; |
599 | > | numAtom = n_atoms; |
600 | #endif | |
601 | < | |
602 | < | for(int i = 0; i < ngroup; i ++){ |
603 | < | mol = &(info->molecules[i]); |
623 | < | numAtom = mol->getNAtoms(); |
624 | < | myAtoms = mol->getMyAtoms(); |
625 | < | mtot = 0.0; |
601 | > | for (int i = 0; i < numAtom; i++) |
602 | > | FglobalGroupMembership.push_back(globalGroupMembership[i] + 1); |
603 | > | |
604 | ||
605 | < | for(int j=0; j < numAtom; j++) |
606 | < | mtot += myAtoms[j]->getMass(); |
607 | < | |
608 | < | for(int j=0; j < numAtom; j++){ |
609 | < | |
610 | < | // We want the local Index: |
611 | < | groupList.push_back(myAtoms[j]->getIndex() + 1); |
634 | < | mfact.push_back(myAtoms[j]->getMass() / mtot); |
605 | > | myMols = info->molecules; |
606 | > | numMol = info->n_mol; |
607 | > | for(int i = 0; i < numMol; i++){ |
608 | > | numCutoffGroups = myMols[i].getNCutoffGroups(); |
609 | > | for(myCutoffGroup =myMols[i].beginCutoffGroup(iterCutoff); |
610 | > | myCutoffGroup != NULL; |
611 | > | myCutoffGroup =myMols[i].nextCutoffGroup(iterCutoff)){ |
612 | ||
613 | < | } |
613 | > | totalMass = myCutoffGroup->getMass(); |
614 | ||
615 | < | groupStart.push_back(curIndex); |
616 | < | curIndex += numAtom; |
617 | < | |
618 | < | } //end for(int i =0 ; i < ngroup; i++) |
615 | > | for(cutoffAtom = myCutoffGroup->beginAtom(iterAtom); |
616 | > | cutoffAtom != NULL; |
617 | > | cutoffAtom = myCutoffGroup->nextAtom(iterAtom)){ |
618 | > | mfact.push_back(cutoffAtom->getMass()/totalMass); |
619 | > | } |
620 | > | } |
621 | } | |
622 | < | else{ |
644 | < | //using atomic cutoff, every single atom is just a group |
645 | < | |
646 | < | #ifdef IS_MPI |
647 | < | ngroup = mpiSim->getMyNlocal(); |
648 | < | #else |
649 | < | ngroup = info->n_atoms; |
650 | < | #endif |
651 | < | |
652 | < | for(int i =0 ; i < ngroup; i++){ |
653 | < | groupStart.push_back(curIndex++); |
654 | < | groupList.push_back((info->atoms[i])->getIndex() + 1); |
655 | < | mfact.push_back(1.0); |
656 | < | |
657 | < | }//end for(int i =0 ; i < ngroup; i++) |
658 | < | |
659 | < | }//end if (info->useMolecularCutoffs) |
660 | < | |
622 | > | |
623 | } |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |