# | Line 14 | Line 14 | using namespace std; | |
---|---|---|
14 | ||
15 | #include "MatVec3.h" | |
16 | ||
17 | + | #include "ConstraintManager.hpp" |
18 | + | |
19 | #ifdef IS_MPI | |
20 | #include "mpiSimulation.hpp" | |
21 | #endif | |
# | Line 42 | Line 44 | SimInfo::SimInfo(){ | |
44 | thermalTime = 0.0; | |
45 | currentTime = 0.0; | |
46 | rCut = 0.0; | |
47 | < | ecr = 0.0; |
46 | < | est = 0.0; |
47 | > | rSw = 0.0; |
48 | ||
49 | haveRcut = 0; | |
50 | < | haveEcr = 0; |
50 | > | haveRsw = 0; |
51 | boxIsInit = 0; | |
52 | ||
53 | resetTime = 1e99; | |
# | Line 63 | Line 64 | SimInfo::SimInfo(){ | |
64 | useReactionField = 0; | |
65 | useGB = 0; | |
66 | useEAM = 0; | |
67 | < | useMolecularCutoffs = 0; |
67 | > | useSolidThermInt = 0; |
68 | > | useLiquidThermInt = 0; |
69 | ||
70 | + | haveCutoffGroups = false; |
71 | + | |
72 | excludes = Exclude::Instance(); | |
73 | ||
74 | myConfiguration = new SimState(); | |
# | Line 74 | Line 78 | SimInfo::SimInfo(){ | |
78 | ||
79 | ngroup = 0; | |
80 | ||
81 | + | consMan = NULL; |
82 | + | |
83 | wrapMeSimInfo( this ); | |
84 | } | |
85 | ||
# | Line 86 | Line 92 | SimInfo::~SimInfo(){ | |
92 | ||
93 | for(i = properties.begin(); i != properties.end(); i++) | |
94 | delete (*i).second; | |
95 | < | |
95 | > | |
96 | > | if (!consMan) |
97 | > | delete consMan; |
98 | } | |
99 | ||
100 | void SimInfo::setBox(double newBox[3]) { | |
# | Line 189 | Line 197 | void SimInfo::calcHmatInv( void ) { | |
197 | ||
198 | if( oldOrtho != orthoRhombic ){ | |
199 | ||
200 | < | if( orthoRhombic ){ |
200 | > | if( orthoRhombic ) { |
201 | sprintf( painCave.errMsg, | |
202 | "OOPSE is switching from the default Non-Orthorhombic\n" | |
203 | "\tto the faster Orthorhombic periodic boundary computations.\n" | |
# | Line 197 | Line 205 | void SimInfo::calcHmatInv( void ) { | |
205 | "\tNon-Orthorhombic computations, make the orthoBoxTolerance\n" | |
206 | "\tvariable ( currently set to %G ) smaller.\n", | |
207 | orthoTolerance); | |
208 | + | painCave.severity = OOPSE_INFO; |
209 | simError(); | |
210 | } | |
211 | else { | |
# | Line 208 | Line 217 | void SimInfo::calcHmatInv( void ) { | |
217 | "\tthe Orthorhombic computations, make the orthoBoxTolerance\n" | |
218 | "\tvariable ( currently set to %G ) larger.\n", | |
219 | orthoTolerance); | |
220 | + | painCave.severity = OOPSE_WARNING; |
221 | simError(); | |
222 | } | |
223 | } | |
# | Line 444 | Line 454 | void SimInfo::refreshSim(){ | |
454 | excl = excludes->getFortranArray(); | |
455 | ||
456 | #ifdef IS_MPI | |
457 | < | n_global = mpiSim->getTotAtoms(); |
457 | > | n_global = mpiSim->getNAtomsGlobal(); |
458 | #else | |
459 | n_global = n_atoms; | |
460 | #endif | |
461 | < | |
461 | > | |
462 | isError = 0; | |
463 | < | |
464 | < | getFortranGroupArray(this, mfact, ngroup, groupList, groupStart); |
465 | < | |
463 | > | |
464 | > | getFortranGroupArrays(this, FglobalGroupMembership, mfact); |
465 | > | //it may not be a good idea to pass the address of first element in vector |
466 | > | //since c++ standard does not require vector to be stored continuously in meomory |
467 | > | //Most of the compilers will organize the memory of vector continuously |
468 | setFsimulation( &fInfo, &n_global, &n_atoms, identArray, &n_exclude, excl, | |
469 | < | &nGlobalExcludes, globalExcludes, molMembershipArray, |
470 | < | &mfact[0], &ngroup, &groupList[0], &groupStart[0], &isError ); |
469 | > | &nGlobalExcludes, globalExcludes, molMembershipArray, |
470 | > | &mfact[0], &ngroup, &FglobalGroupMembership[0], &isError); |
471 | ||
472 | if( isError ){ | |
473 | < | |
473 | > | |
474 | sprintf( painCave.errMsg, | |
475 | < | "There was an error setting the simulation information in fortran.\n" ); |
475 | > | "There was an error setting the simulation information in fortran.\n" ); |
476 | painCave.isFatal = 1; | |
477 | + | painCave.severity = OOPSE_ERROR; |
478 | simError(); | |
479 | } | |
480 | < | |
480 | > | |
481 | #ifdef IS_MPI | |
482 | sprintf( checkPointMsg, | |
483 | "succesfully sent the simulation information to fortran.\n"); | |
484 | MPIcheckPoint(); | |
485 | #endif // is_mpi | |
486 | < | |
486 | > | |
487 | this->ndf = this->getNDF(); | |
488 | this->ndfRaw = this->getNDFraw(); | |
489 | this->ndfTrans = this->getNDFtranslational(); | |
490 | } | |
491 | ||
492 | void SimInfo::setDefaultRcut( double theRcut ){ | |
493 | < | |
493 | > | |
494 | haveRcut = 1; | |
495 | rCut = theRcut; | |
496 | < | |
497 | < | ( rCut > ecr )? rList = rCut + 1.0: rList = ecr + 1.0; |
498 | < | |
486 | < | notifyFortranCutOffs( &rCut, &rList, &ecr, &est ); |
496 | > | rList = rCut + 1.0; |
497 | > | |
498 | > | notifyFortranCutOffs( &rCut, &rSw, &rList ); |
499 | } | |
500 | ||
501 | < | void SimInfo::setDefaultEcr( double theEcr ){ |
501 | > | void SimInfo::setDefaultRcut( double theRcut, double theRsw ){ |
502 | ||
503 | < | haveEcr = 1; |
504 | < | ecr = theEcr; |
493 | < | |
494 | < | ( rCut > ecr )? rList = rCut + 1.0: rList = ecr + 1.0; |
495 | < | |
496 | < | notifyFortranCutOffs( &rCut, &rList, &ecr, &est ); |
503 | > | rSw = theRsw; |
504 | > | setDefaultRcut( theRcut ); |
505 | } | |
506 | ||
499 | – | void SimInfo::setDefaultEcr( double theEcr, double theEst ){ |
507 | ||
501 | – | est = theEst; |
502 | – | setDefaultEcr( theEcr ); |
503 | – | } |
504 | – | |
505 | – | |
508 | void SimInfo::checkCutOffs( void ){ | |
509 | ||
510 | if( boxIsInit ){ | |
# | Line 511 | Line 513 | void SimInfo::checkCutOffs( void ){ | |
513 | ||
514 | if( rCut > maxCutoff ){ | |
515 | sprintf( painCave.errMsg, | |
516 | < | "LJrcut is too large for the current periodic box.\n" |
517 | < | "\tCurrent Value of LJrcut = %G at time %G\n " |
516 | > | "cutoffRadius is too large for the current periodic box.\n" |
517 | > | "\tCurrent Value of cutoffRadius = %G at time %G\n " |
518 | "\tThis is larger than half of at least one of the\n" | |
519 | "\tperiodic box vectors. Right now, the Box matrix is:\n" | |
520 | "\n" | |
# | Line 523 | Line 525 | void SimInfo::checkCutOffs( void ){ | |
525 | Hmat[0][0], Hmat[0][1], Hmat[0][2], | |
526 | Hmat[1][0], Hmat[1][1], Hmat[1][2], | |
527 | Hmat[2][0], Hmat[2][1], Hmat[2][2]); | |
528 | + | painCave.severity = OOPSE_ERROR; |
529 | painCave.isFatal = 1; | |
530 | simError(); | |
531 | < | } |
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 | < | } |
531 | > | } |
532 | } else { | |
533 | // initialize this stuff before using it, OK? | |
534 | sprintf( painCave.errMsg, | |
535 | "Trying to check cutoffs without a box.\n" | |
536 | "\tOOPSE should have better programmers than that.\n" ); | |
537 | + | painCave.severity = OOPSE_ERROR; |
538 | painCave.isFatal = 1; | |
539 | simError(); | |
540 | } | |
# | Line 595 | Line 578 | GenericData* SimInfo::getProperty(const string& propNa | |
578 | } | |
579 | ||
580 | ||
581 | < | void getFortranGroupArray(SimInfo* info, vector<double>& mfact, int& ngroup, |
582 | < | vector<int>& groupList, vector<int>& groupStart){ |
583 | < | Molecule* mol; |
581 | > | void SimInfo::getFortranGroupArrays(SimInfo* info, |
582 | > | vector<int>& FglobalGroupMembership, |
583 | > | vector<double>& mfact){ |
584 | > | |
585 | > | Molecule* myMols; |
586 | > | Atom** myAtoms; |
587 | int numAtom; | |
588 | < | int curIndex; |
589 | < | |
588 | > | double mtot; |
589 | > | int numMol; |
590 | > | int numCutoffGroups; |
591 | > | CutoffGroup* myCutoffGroup; |
592 | > | vector<CutoffGroup*>::iterator iterCutoff; |
593 | > | Atom* cutoffAtom; |
594 | > | vector<Atom*>::iterator iterAtom; |
595 | > | int atomIndex; |
596 | > | double totalMass; |
597 | > | |
598 | mfact.clear(); | |
599 | < | groupList.clear(); |
600 | < | groupStart.clear(); |
599 | > | FglobalGroupMembership.clear(); |
600 | > | |
601 | ||
602 | < | //Be careful, fortran array begin at 1 |
609 | < | curIndex = 1; |
610 | < | |
611 | < | if(info->useMolecularCutoffs){ |
612 | < | //if using molecular cutoff |
613 | < | ngroup = info->n_mol; |
614 | < | |
615 | < | for(int i = 0; i < ngroup; i ++){ |
616 | < | mol = &(info->molecules[i]); |
617 | < | numAtom = mol->getNAtoms(); |
618 | < | |
619 | < | for(int j=0; j < numAtom; j++){ |
602 | > | // Fix the silly fortran indexing problem |
603 | #ifdef IS_MPI | |
604 | < | groupList.push_back((info->atoms[i])->getGlobalIndex() + 1); |
604 | > | numAtom = mpiSim->getNAtomsGlobal(); |
605 | #else | |
606 | < | groupList.push_back((info->atoms[i])->getIndex() + 1); |
606 | > | numAtom = n_atoms; |
607 | #endif | |
608 | < | }//for(int j=0; j < numAtom; j++) |
609 | < | |
610 | < | groupStart.push_back(curIndex); |
611 | < | curIndex += numAtom; |
608 | > | for (int i = 0; i < numAtom; i++) |
609 | > | FglobalGroupMembership.push_back(globalGroupMembership[i] + 1); |
610 | > | |
611 | > | |
612 | > | myMols = info->molecules; |
613 | > | numMol = info->n_mol; |
614 | > | for(int i = 0; i < numMol; i++){ |
615 | > | numCutoffGroups = myMols[i].getNCutoffGroups(); |
616 | > | for(myCutoffGroup =myMols[i].beginCutoffGroup(iterCutoff); |
617 | > | myCutoffGroup != NULL; |
618 | > | myCutoffGroup =myMols[i].nextCutoffGroup(iterCutoff)){ |
619 | > | |
620 | > | totalMass = myCutoffGroup->getMass(); |
621 | ||
622 | < | }//end for(int i =0 ; i < ngroup; i++) |
622 | > | for(cutoffAtom = myCutoffGroup->beginAtom(iterAtom); |
623 | > | cutoffAtom != NULL; |
624 | > | cutoffAtom = myCutoffGroup->nextAtom(iterAtom)){ |
625 | > | mfact.push_back(cutoffAtom->getMass()/totalMass); |
626 | > | } |
627 | > | } |
628 | } | |
632 | – | else{ |
633 | – | //using atomic cutoff, every single atom is just a group |
634 | – | ngroup = info->n_atoms; |
635 | – | for(int i =0 ; i < ngroup; i++){ |
636 | – | groupStart.push_back(curIndex++); |
629 | ||
638 | – | #ifdef IS_MPI |
639 | – | groupList.push_back((info->atoms[i])->getGlobalIndex() + 1); |
640 | – | #else |
641 | – | groupList.push_back((info->atoms[i])->getIndex() + 1); |
642 | – | #endif |
643 | – | |
644 | – | }//end for(int i =0 ; i < ngroup; i++) |
645 | – | |
646 | – | }//end if (info->useMolecularCutoffs) |
647 | – | |
630 | } |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |