| 73 |
|
|
| 74 |
|
void ForceFields::doForces( int calcPot, int calcStress ){ |
| 75 |
|
|
| 76 |
< |
int i, isError; |
| 76 |
> |
int i, j, isError; |
| 77 |
|
double* frc; |
| 78 |
|
double* pos; |
| 79 |
|
double* trq; |
| 80 |
|
double* A; |
| 81 |
< |
double* u_l;; |
| 81 |
> |
double* u_l; |
| 82 |
> |
double* rc; |
| 83 |
> |
double* massRatio; |
| 84 |
|
SimState* config; |
| 85 |
|
|
| 86 |
+ |
Molecule* myMols; |
| 87 |
+ |
Atom** myAtoms; |
| 88 |
+ |
int numAtom; |
| 89 |
+ |
int curIndex; |
| 90 |
+ |
double mtot; |
| 91 |
+ |
int numMol; |
| 92 |
+ |
int numCutoffGroups; |
| 93 |
+ |
CutoffGroup* myCutoffGroup; |
| 94 |
+ |
vector<CutoffGroup*>::iterator iterCutoff; |
| 95 |
+ |
Atom* cutoffAtom; |
| 96 |
+ |
vector<Atom*>::iterator iterAtom; |
| 97 |
+ |
double com[3]; |
| 98 |
+ |
double tempPos[3]; |
| 99 |
+ |
int atomIndex; |
| 100 |
+ |
|
| 101 |
|
short int passedCalcPot = (short int)calcPot; |
| 102 |
|
short int passedCalcStress = (short int)calcStress; |
| 103 |
|
|
| 115 |
|
for(i=0; i<entry_plug->n_mol; i++ ){ |
| 116 |
|
// CalcForces in molecules takes care of mapping rigid body coordinates |
| 117 |
|
// into atomic coordinates |
| 118 |
< |
entry_plug->molecules[i].calcForces(); |
| 118 |
> |
entry_plug->molecules[i].calcForces(); |
| 119 |
|
} |
| 120 |
|
|
| 121 |
|
#ifdef PROFILE |
| 130 |
|
A = config->getAmatArray(); |
| 131 |
|
u_l = config->getUlArray(); |
| 132 |
|
|
| 133 |
+ |
if(entry_plug->haveCutoffGroups){ |
| 134 |
+ |
//if |
| 135 |
+ |
myMols = entry_plug->molecules; |
| 136 |
+ |
numMol = entry_plug->n_mol; |
| 137 |
+ |
for(int i = 0; i < numMol; i++){ |
| 138 |
+ |
numAtom = myMols[i].getNAtoms(); |
| 139 |
+ |
myAtoms = myMols[i].getMyAtoms(); |
| 140 |
+ |
|
| 141 |
+ |
|
| 142 |
+ |
for(int j = 0; j < numAtom; j++){ |
| 143 |
+ |
#ifdef IS_MPI |
| 144 |
+ |
atomIndex = myAtoms[j]->getGlobalIndex(); |
| 145 |
+ |
#else |
| 146 |
+ |
atomIndex = myAtoms[j]->getIndex(); |
| 147 |
+ |
#endif |
| 148 |
+ |
|
| 149 |
+ |
if(myMols[i].belongToCutoffGroup(atomIndex)) |
| 150 |
+ |
continue; |
| 151 |
+ |
else{ |
| 152 |
+ |
myAtoms[j]->getPos(tempPos); |
| 153 |
+ |
myAtoms[j]->setRc(tempPos); |
| 154 |
+ |
} |
| 155 |
+ |
|
| 156 |
+ |
} |
| 157 |
+ |
|
| 158 |
+ |
numCutoffGroups = myMols[i].getNCutoffGroups(); |
| 159 |
+ |
for(myCutoffGroup =myMols[i].beginCutoffGroup(iterCutoff); myCutoffGroup != NULL; |
| 160 |
+ |
myCutoffGroup =myMols[i].nextCutoffGroup(iterCutoff)){ |
| 161 |
+ |
//get center of mass of the cutoff group |
| 162 |
+ |
myCutoffGroup->getCOM(com); |
| 163 |
+ |
|
| 164 |
+ |
for(cutoffAtom = myCutoffGroup->beginAtom(iterAtom); cutoffAtom != NULL; |
| 165 |
+ |
cutoffAtom = myCutoffGroup->beginAtom(iterAtom)){ |
| 166 |
+ |
cutoffAtom->setRc(com); |
| 167 |
+ |
} |
| 168 |
+ |
|
| 169 |
+ |
}// end for(myCutoffGroup) |
| 170 |
+ |
|
| 171 |
+ |
}//end for(int i = 0) |
| 172 |
+ |
|
| 173 |
+ |
rc = config->getRcArray(); |
| 174 |
+ |
} |
| 175 |
+ |
else{ |
| 176 |
+ |
// center of mass of the group is the same as position of the atom if cutoff group does not exist |
| 177 |
+ |
rc = pos; |
| 178 |
+ |
} |
| 179 |
+ |
|
| 180 |
+ |
|
| 181 |
+ |
|
| 182 |
|
isError = 0; |
| 183 |
|
entry_plug->lrPot = 0.0; |
| 184 |
|
|
| 192 |
|
#endif |
| 193 |
|
|
| 194 |
|
fortranForceLoop( pos, |
| 195 |
+ |
rc, |
| 196 |
|
A, |
| 197 |
|
u_l, |
| 198 |
|
frc, |
| 203 |
|
&passedCalcStress, |
| 204 |
|
&isError ); |
| 205 |
|
|
| 206 |
+ |
|
| 207 |
|
#ifdef PROFILE |
| 208 |
|
endProfile(pro8); |
| 209 |
|
#endif |
