# | Line 9 | Line 9 | |
---|---|---|
9 | #include "parse_me.h" | |
10 | #include "Integrator.hpp" | |
11 | #include "simError.h" | |
12 | + | #include "RigidBody.hpp" |
13 | + | //#include "ConjugateMinimizer.hpp" |
14 | + | #include "OOPSEMinimizer.hpp" |
15 | ||
16 | #ifdef IS_MPI | |
17 | #include "mpiBASS.h" | |
# | Line 24 | Line 27 | |
27 | #define NPTxyz_ENS 4 | |
28 | ||
29 | ||
30 | < | #define FF_DUFF 0 |
31 | < | #define FF_LJ 1 |
32 | < | #define FF_EAM 2 |
30 | > | #define FF_DUFF 0 |
31 | > | #define FF_LJ 1 |
32 | > | #define FF_EAM 2 |
33 | > | #define FF_H2O 3 |
34 | ||
35 | using namespace std; | |
36 | ||
# | Line 144 | Line 148 | void SimSetup::createSim(void){ | |
148 | ||
149 | makeOutNames(); | |
150 | ||
151 | < | // make the integrator |
152 | < | |
153 | < | makeIntegrator(); |
154 | < | |
151 | > | if (globals->haveMinimizer()) |
152 | > | // make minimizer |
153 | > | makeMinimizer(); |
154 | > | else |
155 | > | // make the integrator |
156 | > | makeIntegrator(); |
157 | > | |
158 | #ifdef IS_MPI | |
159 | mpiSim->mpiRefresh(); | |
160 | #endif | |
# | Line 159 | Line 166 | void SimSetup::makeMolecules(void){ | |
166 | ||
167 | ||
168 | void SimSetup::makeMolecules(void){ | |
169 | < | int k; |
170 | < | int i, j, exI, exJ, tempEx, stampID, atomOffset, excludeOffset; |
169 | > | int i, j, k; |
170 | > | int exI, exJ, exK, exL, slI; |
171 | > | int tempI, tempJ, tempK, tempL; |
172 | > | int molI; |
173 | > | int stampID, atomOffset, rbOffset; |
174 | molInit molInfo; | |
175 | DirectionalAtom* dAtom; | |
176 | + | RigidBody* myRB; |
177 | + | StuntDouble* mySD; |
178 | LinkedAssign* extras; | |
179 | LinkedAssign* current_extra; | |
180 | AtomStamp* currentAtom; | |
181 | BondStamp* currentBond; | |
182 | BendStamp* currentBend; | |
183 | TorsionStamp* currentTorsion; | |
184 | + | RigidBodyStamp* currentRigidBody; |
185 | ||
186 | bond_pair* theBonds; | |
187 | bend_set* theBends; | |
188 | torsion_set* theTorsions; | |
189 | ||
190 | + | set<int> skipList; |
191 | ||
192 | + | double phi, theta, psi; |
193 | + | |
194 | //init the forceField paramters | |
195 | ||
196 | the_ff->readParams(); | |
197 | ||
182 | – | |
198 | // init the atoms | |
199 | ||
200 | < | double ux, uy, uz, u, uSqr; |
200 | > | int nMembers, nNew, rb1, rb2; |
201 | ||
202 | for (k = 0; k < nInfo; k++){ | |
203 | the_ff->setSimInfo(&(info[k])); | |
204 | ||
205 | atomOffset = 0; | |
206 | < | excludeOffset = 0; |
206 | > | |
207 | for (i = 0; i < info[k].n_mol; i++){ | |
208 | stampID = info[k].molecules[i].getStampID(); | |
209 | ||
# | Line 196 | Line 211 | void SimSetup::makeMolecules(void){ | |
211 | molInfo.nBonds = comp_stamps[stampID]->getNBonds(); | |
212 | molInfo.nBends = comp_stamps[stampID]->getNBends(); | |
213 | molInfo.nTorsions = comp_stamps[stampID]->getNTorsions(); | |
214 | < | molInfo.nExcludes = molInfo.nBonds + molInfo.nBends + molInfo.nTorsions; |
215 | < | |
214 | > | molInfo.nRigidBodies = comp_stamps[stampID]->getNRigidBodies(); |
215 | > | |
216 | molInfo.myAtoms = &(info[k].atoms[atomOffset]); | |
202 | – | molInfo.myExcludes = &(info[k].excludes[excludeOffset]); |
203 | – | molInfo.myBonds = new Bond * [molInfo.nBonds]; |
204 | – | molInfo.myBends = new Bend * [molInfo.nBends]; |
205 | – | molInfo.myTorsions = new Torsion * [molInfo.nTorsions]; |
217 | ||
218 | + | if (molInfo.nBonds > 0) |
219 | + | molInfo.myBonds = new (Bond *) [molInfo.nBonds]; |
220 | + | else |
221 | + | molInfo.myBonds = NULL; |
222 | + | |
223 | + | if (molInfo.nBends > 0) |
224 | + | molInfo.myBends = new (Bend *) [molInfo.nBends]; |
225 | + | else |
226 | + | molInfo.myBends = NULL; |
227 | + | |
228 | + | if (molInfo.nTorsions > 0) |
229 | + | molInfo.myTorsions = new (Torsion *) [molInfo.nTorsions]; |
230 | + | else |
231 | + | molInfo.myTorsions = NULL; |
232 | + | |
233 | theBonds = new bond_pair[molInfo.nBonds]; | |
234 | theBends = new bend_set[molInfo.nBends]; | |
235 | theTorsions = new torsion_set[molInfo.nTorsions]; | |
236 | < | |
236 | > | |
237 | // make the Atoms | |
238 | ||
239 | for (j = 0; j < molInfo.nAtoms; j++){ | |
240 | currentAtom = comp_stamps[stampID]->getAtom(j); | |
241 | + | |
242 | if (currentAtom->haveOrientation()){ | |
243 | dAtom = new DirectionalAtom((j + atomOffset), | |
244 | info[k].getConfiguration()); | |
245 | info[k].n_oriented++; | |
246 | molInfo.myAtoms[j] = dAtom; | |
247 | ||
248 | < | ux = currentAtom->getOrntX(); |
249 | < | uy = currentAtom->getOrntY(); |
250 | < | uz = currentAtom->getOrntZ(); |
248 | > | // Directional Atoms have standard unit vectors which are oriented |
249 | > | // in space using the three Euler angles. We assume the standard |
250 | > | // unit vector was originally along the z axis below. |
251 | ||
252 | < | uSqr = (ux * ux) + (uy * uy) + (uz * uz); |
252 | > | phi = currentAtom->getEulerPhi() * M_PI / 180.0; |
253 | > | theta = currentAtom->getEulerTheta() * M_PI / 180.0; |
254 | > | psi = currentAtom->getEulerPsi()* M_PI / 180.0; |
255 | ||
256 | < | u = sqrt(uSqr); |
257 | < | ux = ux / u; |
229 | < | uy = uy / u; |
230 | < | uz = uz / u; |
231 | < | |
232 | < | dAtom->setSUx(ux); |
233 | < | dAtom->setSUy(uy); |
234 | < | dAtom->setSUz(uz); |
256 | > | dAtom->setUnitFrameFromEuler(phi, theta, psi); |
257 | > | |
258 | } | |
259 | else{ | |
260 | < | molInfo.myAtoms[j] = new GeneralAtom((j + atomOffset), |
261 | < | info[k].getConfiguration()); |
260 | > | |
261 | > | molInfo.myAtoms[j] = new Atom((j + atomOffset), info[k].getConfiguration()); |
262 | } | |
263 | + | |
264 | molInfo.myAtoms[j]->setType(currentAtom->getType()); | |
265 | ||
266 | #ifdef IS_MPI | |
# | Line 252 | Line 276 | void SimSetup::makeMolecules(void){ | |
276 | theBonds[j].a = currentBond->getA() + atomOffset; | |
277 | theBonds[j].b = currentBond->getB() + atomOffset; | |
278 | ||
279 | < | exI = theBonds[j].a; |
280 | < | exJ = theBonds[j].b; |
279 | > | tempI = theBonds[j].a; |
280 | > | tempJ = theBonds[j].b; |
281 | ||
258 | – | // exclude_I must always be the smaller of the pair |
259 | – | if (exI > exJ){ |
260 | – | tempEx = exI; |
261 | – | exI = exJ; |
262 | – | exJ = tempEx; |
263 | – | } |
282 | #ifdef IS_MPI | |
283 | < | tempEx = exI; |
284 | < | exI = info[k].atoms[tempEx]->getGlobalIndex() + 1; |
285 | < | tempEx = exJ; |
286 | < | exJ = info[k].atoms[tempEx]->getGlobalIndex() + 1; |
283 | > | exI = info[k].atoms[tempI]->getGlobalIndex() + 1; |
284 | > | exJ = info[k].atoms[tempJ]->getGlobalIndex() + 1; |
285 | > | #else |
286 | > | exI = tempI + 1; |
287 | > | exJ = tempJ + 1; |
288 | > | #endif |
289 | ||
290 | < | info[k].excludes[j + excludeOffset]->setPair(exI, exJ); |
271 | < | #else // isn't MPI |
272 | < | |
273 | < | info[k].excludes[j + excludeOffset]->setPair((exI + 1), (exJ + 1)); |
274 | < | #endif //is_mpi |
290 | > | info[k].excludes->addPair(exI, exJ); |
291 | } | |
276 | – | excludeOffset += molInfo.nBonds; |
292 | ||
293 | //make the bends | |
294 | for (j = 0; j < molInfo.nBends; j++){ | |
# | Line 321 | Line 336 | void SimSetup::makeMolecules(void){ | |
336 | ||
337 | current_extra = current_extra->getNext(); | |
338 | } | |
324 | – | } |
325 | – | |
326 | – | if (!theBends[j].isGhost){ |
327 | – | exI = theBends[j].a; |
328 | – | exJ = theBends[j].c; |
339 | } | |
330 | – | else{ |
331 | – | exI = theBends[j].a; |
332 | – | exJ = theBends[j].b; |
333 | – | } |
340 | ||
341 | < | // exclude_I must always be the smaller of the pair |
342 | < | if (exI > exJ){ |
343 | < | tempEx = exI; |
344 | < | exI = exJ; |
345 | < | exJ = tempEx; |
340 | < | } |
341 | > | if (theBends[j].isGhost) { |
342 | > | |
343 | > | tempI = theBends[j].a; |
344 | > | tempJ = theBends[j].b; |
345 | > | |
346 | #ifdef IS_MPI | |
347 | < | tempEx = exI; |
348 | < | exI = info[k].atoms[tempEx]->getGlobalIndex() + 1; |
349 | < | tempEx = exJ; |
350 | < | exJ = info[k].atoms[tempEx]->getGlobalIndex() + 1; |
347 | > | exI = info[k].atoms[tempI]->getGlobalIndex() + 1; |
348 | > | exJ = info[k].atoms[tempJ]->getGlobalIndex() + 1; |
349 | > | #else |
350 | > | exI = tempI + 1; |
351 | > | exJ = tempJ + 1; |
352 | > | #endif |
353 | > | info[k].excludes->addPair(exI, exJ); |
354 | ||
355 | < | info[k].excludes[j + excludeOffset]->setPair(exI, exJ); |
356 | < | #else // isn't MPI |
357 | < | info[k].excludes[j + excludeOffset]->setPair((exI + 1), (exJ + 1)); |
358 | < | #endif //is_mpi |
355 | > | } else { |
356 | > | |
357 | > | tempI = theBends[j].a; |
358 | > | tempJ = theBends[j].b; |
359 | > | tempK = theBends[j].c; |
360 | > | |
361 | > | #ifdef IS_MPI |
362 | > | exI = info[k].atoms[tempI]->getGlobalIndex() + 1; |
363 | > | exJ = info[k].atoms[tempJ]->getGlobalIndex() + 1; |
364 | > | exK = info[k].atoms[tempK]->getGlobalIndex() + 1; |
365 | > | #else |
366 | > | exI = tempI + 1; |
367 | > | exJ = tempJ + 1; |
368 | > | exK = tempK + 1; |
369 | > | #endif |
370 | > | |
371 | > | info[k].excludes->addPair(exI, exK); |
372 | > | info[k].excludes->addPair(exI, exJ); |
373 | > | info[k].excludes->addPair(exJ, exK); |
374 | > | } |
375 | } | |
352 | – | excludeOffset += molInfo.nBends; |
376 | ||
377 | for (j = 0; j < molInfo.nTorsions; j++){ | |
378 | currentTorsion = comp_stamps[stampID]->getTorsion(j); | |
# | Line 358 | Line 381 | void SimSetup::makeMolecules(void){ | |
381 | theTorsions[j].c = currentTorsion->getC() + atomOffset; | |
382 | theTorsions[j].d = currentTorsion->getD() + atomOffset; | |
383 | ||
384 | < | exI = theTorsions[j].a; |
385 | < | exJ = theTorsions[j].d; |
384 | > | tempI = theTorsions[j].a; |
385 | > | tempJ = theTorsions[j].b; |
386 | > | tempK = theTorsions[j].c; |
387 | > | tempL = theTorsions[j].d; |
388 | ||
364 | – | // exclude_I must always be the smaller of the pair |
365 | – | if (exI > exJ){ |
366 | – | tempEx = exI; |
367 | – | exI = exJ; |
368 | – | exJ = tempEx; |
369 | – | } |
389 | #ifdef IS_MPI | |
390 | < | tempEx = exI; |
391 | < | exI = info[k].atoms[tempEx]->getGlobalIndex() + 1; |
392 | < | tempEx = exJ; |
393 | < | exJ = info[k].atoms[tempEx]->getGlobalIndex() + 1; |
390 | > | exI = info[k].atoms[tempI]->getGlobalIndex() + 1; |
391 | > | exJ = info[k].atoms[tempJ]->getGlobalIndex() + 1; |
392 | > | exK = info[k].atoms[tempK]->getGlobalIndex() + 1; |
393 | > | exL = info[k].atoms[tempL]->getGlobalIndex() + 1; |
394 | > | #else |
395 | > | exI = tempI + 1; |
396 | > | exJ = tempJ + 1; |
397 | > | exK = tempK + 1; |
398 | > | exL = tempL + 1; |
399 | > | #endif |
400 | ||
401 | < | info[k].excludes[j + excludeOffset]->setPair(exI, exJ); |
402 | < | #else // isn't MPI |
403 | < | info[k].excludes[j + excludeOffset]->setPair((exI + 1), (exJ + 1)); |
404 | < | #endif //is_mpi |
401 | > | info[k].excludes->addPair(exI, exJ); |
402 | > | info[k].excludes->addPair(exI, exK); |
403 | > | info[k].excludes->addPair(exI, exL); |
404 | > | info[k].excludes->addPair(exJ, exK); |
405 | > | info[k].excludes->addPair(exJ, exL); |
406 | > | info[k].excludes->addPair(exK, exL); |
407 | } | |
381 | – | excludeOffset += molInfo.nTorsions; |
408 | ||
409 | + | for (j = 0; j < molInfo.nRigidBodies; j++){ |
410 | ||
411 | < | // send the arrays off to the forceField for init. |
411 | > | currentRigidBody = comp_stamps[stampID]->getRigidBody(j); |
412 | > | nMembers = currentRigidBody->getNMembers(); |
413 | > | |
414 | > | // Create the Rigid Body: |
415 | > | |
416 | > | myRB = new RigidBody(); |
417 | > | |
418 | > | for (rb1 = 0; rb1 < nMembers; rb1++) { |
419 | > | |
420 | > | // molI is atom numbering inside this molecule |
421 | > | molI = currentRigidBody->getMember(rb1); |
422 | > | |
423 | > | // tempI is atom numbering on local processor |
424 | > | tempI = molI + atomOffset; |
425 | > | |
426 | > | // currentAtom is the AtomStamp (which we need for |
427 | > | // rigid body reference positions) |
428 | > | currentAtom = comp_stamps[stampID]->getAtom(molI); |
429 | ||
430 | + | // When we add to the rigid body, add the atom itself and |
431 | + | // the stamp info: |
432 | + | |
433 | + | myRB->addAtom(info[k].atoms[tempI], currentAtom); |
434 | + | |
435 | + | // Add this atom to the Skip List for the integrators |
436 | + | #ifdef IS_MPI |
437 | + | slI = info[k].atoms[tempI]->getGlobalIndex(); |
438 | + | #else |
439 | + | slI = tempI; |
440 | + | #endif |
441 | + | skipList.insert(slI); |
442 | + | |
443 | + | } |
444 | + | |
445 | + | for(rb1 = 0; rb1 < nMembers - 1; rb1++) { |
446 | + | for(rb2 = rb1+1; rb2 < nMembers; rb2++) { |
447 | + | |
448 | + | tempI = currentRigidBody->getMember(rb1); |
449 | + | tempJ = currentRigidBody->getMember(rb2); |
450 | + | |
451 | + | // Some explanation is required here. |
452 | + | // Fortran indexing starts at 1, while c indexing starts at 0 |
453 | + | // Also, in parallel computations, the GlobalIndex is |
454 | + | // used for the exclude list: |
455 | + | |
456 | + | #ifdef IS_MPI |
457 | + | exI = info[k].atoms[tempI]->getGlobalIndex() + 1; |
458 | + | exJ = info[k].atoms[tempJ]->getGlobalIndex() + 1; |
459 | + | #else |
460 | + | exI = tempI + 1; |
461 | + | exJ = tempJ + 1; |
462 | + | #endif |
463 | + | |
464 | + | info[k].excludes->addPair(exI, exJ); |
465 | + | |
466 | + | } |
467 | + | } |
468 | + | } |
469 | + | |
470 | + | // send the arrays off to the forceField for init. |
471 | + | |
472 | the_ff->initializeAtoms(molInfo.nAtoms, molInfo.myAtoms); | |
473 | the_ff->initializeBonds(molInfo.nBonds, molInfo.myBonds, theBonds); | |
474 | the_ff->initializeBends(molInfo.nBends, molInfo.myBends, theBends); | |
475 | the_ff->initializeTorsions(molInfo.nTorsions, molInfo.myTorsions, | |
476 | theTorsions); | |
477 | ||
392 | – | |
478 | info[k].molecules[i].initialize(molInfo); | |
479 | ||
480 | ||
# | Line 398 | Line 483 | void SimSetup::makeMolecules(void){ | |
483 | delete[] theBends; | |
484 | delete[] theTorsions; | |
485 | } | |
486 | + | |
487 | + | // build up the integrableObjects vector: |
488 | + | |
489 | + | for (i = 0; i < info[k].n_atoms; i++) { |
490 | + | |
491 | + | #ifdef IS_MPI |
492 | + | slI = info[k].atoms[i]->getGlobalIndex(); |
493 | + | #else |
494 | + | slI = i; |
495 | + | #endif |
496 | + | |
497 | + | if (skipList.find(slI) == skipList.end()) { |
498 | + | mySD = (StuntDouble *) info[k].atoms[i]; |
499 | + | info[k].integrableObjects.push_back(mySD); |
500 | + | } |
501 | + | } |
502 | + | for (i = 0; i < info[k].rigidBodies.size(); i++) { |
503 | + | mySD = (StuntDouble *) info[k].rigidBodies[i]; |
504 | + | info[k].integrableObjects.push_back(mySD); |
505 | + | } |
506 | + | |
507 | } | |
508 | ||
509 | #ifdef IS_MPI | |
# | Line 407 | Line 513 | void SimSetup::makeMolecules(void){ | |
513 | ||
514 | // clean up the forcefield | |
515 | ||
516 | < | the_ff->calcRcut(); |
516 | > | if (!globals->haveLJrcut()){ |
517 | > | |
518 | > | the_ff->calcRcut(); |
519 | > | |
520 | > | } else { |
521 | > | |
522 | > | the_ff->setRcut( globals->getLJrcut() ); |
523 | > | } |
524 | > | |
525 | the_ff->cleanMe(); | |
526 | } | |
527 | ||
# | Line 608 | Line 722 | void SimSetup::gatherInfo(void){ | |
722 | } | |
723 | else if (!strcasecmp(force_field, "EAM")){ | |
724 | ffCase = FF_EAM; | |
725 | + | } |
726 | + | else if (!strcasecmp(force_field, "WATER")){ |
727 | + | ffCase = FF_H2O; |
728 | } | |
729 | else{ | |
730 | sprintf(painCave.errMsg, "SimSetup Error. Unrecognized force field -> %s\n", | |
# | Line 811 | Line 928 | void SimSetup::gatherInfo(void){ | |
928 | for (int i = 0; i < nInfo; i++){ | |
929 | info[i].setSeed(seedValue); | |
930 | } | |
931 | < | |
931 | > | |
932 | #ifdef IS_MPI | |
933 | < | strcpy(checkPointMsg, "Succesfully gathered all information from Bass\n"); |
933 | > | strcpy(checkPointMsg, "Successfully gathered all information from Bass\n"); |
934 | MPIcheckPoint(); | |
935 | #endif // is_mpi | |
936 | } | |
# | Line 1110 | Line 1227 | void SimSetup::createFF(void){ | |
1227 | the_ff = new EAM_FF(); | |
1228 | break; | |
1229 | ||
1230 | + | case FF_H2O: |
1231 | + | the_ff = new WATER(); |
1232 | + | break; |
1233 | + | |
1234 | default: | |
1235 | sprintf(painCave.errMsg, | |
1236 | "SimSetup Error. Unrecognized force field in case statement.\n"); | |
# | Line 1182 | Line 1303 | void SimSetup::calcSysValues(void){ | |
1303 | tot_bonds = 0; | |
1304 | tot_bends = 0; | |
1305 | tot_torsions = 0; | |
1306 | + | tot_rigid = 0; |
1307 | for (i = 0; i < n_components; i++){ | |
1308 | tot_atoms += components_nmol[i] * comp_stamps[i]->getNAtoms(); | |
1309 | tot_bonds += components_nmol[i] * comp_stamps[i]->getNBonds(); | |
1310 | tot_bends += components_nmol[i] * comp_stamps[i]->getNBends(); | |
1311 | tot_torsions += components_nmol[i] * comp_stamps[i]->getNTorsions(); | |
1312 | + | tot_rigid += components_nmol[i] * comp_stamps[i]->getNRigidBodies(); |
1313 | } | |
1314 | < | |
1314 | > | |
1315 | tot_SRI = tot_bonds + tot_bends + tot_torsions; | |
1316 | molMembershipArray = new int[tot_atoms]; | |
1317 | ||
# | Line 1210 | Line 1333 | void SimSetup::mpiMolDivide(void){ | |
1333 | int i, j, k; | |
1334 | int localMol, allMol; | |
1335 | int local_atoms, local_bonds, local_bends, local_torsions, local_SRI; | |
1336 | + | int local_rigid; |
1337 | ||
1338 | mpiSim = new mpiSimulation(info); | |
1339 | ||
# | Line 1226 | Line 1350 | void SimSetup::mpiMolDivide(void){ | |
1350 | local_bonds = 0; | |
1351 | local_bends = 0; | |
1352 | local_torsions = 0; | |
1353 | + | local_rigid = 0; |
1354 | globalAtomIndex = 0; | |
1355 | ||
1231 | – | |
1356 | for (i = 0; i < n_components; i++){ | |
1357 | for (j = 0; j < components_nmol[i]; j++){ | |
1358 | if (mol2proc[allMol] == worldRank){ | |
# | Line 1236 | Line 1360 | void SimSetup::mpiMolDivide(void){ | |
1360 | local_bonds += comp_stamps[i]->getNBonds(); | |
1361 | local_bends += comp_stamps[i]->getNBends(); | |
1362 | local_torsions += comp_stamps[i]->getNTorsions(); | |
1363 | + | local_rigid += comp_stamps[i]->getNRigidBodies(); |
1364 | localMol++; | |
1365 | } | |
1366 | for (k = 0; k < comp_stamps[i]->getNAtoms(); k++){ | |
# | Line 1249 | Line 1374 | void SimSetup::mpiMolDivide(void){ | |
1374 | local_SRI = local_bonds + local_bends + local_torsions; | |
1375 | ||
1376 | info[0].n_atoms = mpiSim->getMyNlocal(); | |
1377 | + | |
1378 | ||
1379 | if (local_atoms != info[0].n_atoms){ | |
1380 | sprintf(painCave.errMsg, | |
# | Line 1281 | Line 1407 | void SimSetup::makeSysArrays(void){ | |
1407 | ||
1408 | Atom** the_atoms; | |
1409 | Molecule* the_molecules; | |
1284 | – | Exclude** the_excludes; |
1410 | ||
1286 | – | |
1411 | for (l = 0; l < nInfo; l++){ | |
1412 | // create the atom and short range interaction arrays | |
1413 | ||
# | Line 1326 | Line 1450 | void SimSetup::makeSysArrays(void){ | |
1450 | ||
1451 | #endif // is_mpi | |
1452 | ||
1453 | < | |
1454 | < | if (info[l].n_SRI){ |
1455 | < | Exclude::createArray(info[l].n_SRI); |
1332 | < | the_excludes = new Exclude * [info[l].n_SRI]; |
1333 | < | for (int ex = 0; ex < info[l].n_SRI; ex++){ |
1334 | < | the_excludes[ex] = new Exclude(ex); |
1335 | < | } |
1336 | < | info[l].globalExcludes = new int; |
1337 | < | info[l].n_exclude = info[l].n_SRI; |
1338 | < | } |
1339 | < | else{ |
1340 | < | Exclude::createArray(1); |
1341 | < | the_excludes = new Exclude * ; |
1342 | < | the_excludes[0] = new Exclude(0); |
1343 | < | the_excludes[0]->setPair(0, 0); |
1344 | < | info[l].globalExcludes = new int; |
1345 | < | info[l].globalExcludes[0] = 0; |
1346 | < | info[l].n_exclude = 0; |
1347 | < | } |
1348 | < | |
1453 | > | info[l].globalExcludes = new int; |
1454 | > | info[l].globalExcludes[0] = 0; |
1455 | > | |
1456 | // set the arrays into the SimInfo object | |
1457 | ||
1458 | info[l].atoms = the_atoms; | |
1459 | info[l].molecules = the_molecules; | |
1460 | info[l].nGlobalExcludes = 0; | |
1354 | – | info[l].excludes = the_excludes; |
1461 | ||
1462 | the_ff->setSimInfo(info); | |
1463 | } | |
# | Line 1625 | Line 1731 | void SimSetup::setupZConstraint(SimInfo& theInfo){ | |
1731 | ||
1732 | theInfo.addProperty(zconsForcePolicy); | |
1733 | ||
1734 | + | //set zcons gap |
1735 | + | DoubleData* zconsGap = new DoubleData(); |
1736 | + | zconsGap->setID(ZCONSGAP_ID); |
1737 | + | |
1738 | + | if (globals->haveZConsGap()){ |
1739 | + | zconsGap->setData(globals->getZconsGap()); |
1740 | + | theInfo.addProperty(zconsGap); |
1741 | + | } |
1742 | + | |
1743 | + | //set zcons fixtime |
1744 | + | DoubleData* zconsFixtime = new DoubleData(); |
1745 | + | zconsFixtime->setID(ZCONSFIXTIME_ID); |
1746 | + | |
1747 | + | if (globals->haveZConsFixTime()){ |
1748 | + | zconsFixtime->setData(globals->getZconsFixtime()); |
1749 | + | theInfo.addProperty(zconsFixtime); |
1750 | + | } |
1751 | + | |
1752 | + | //set zconsUsingSMD |
1753 | + | IntData* zconsUsingSMD = new IntData(); |
1754 | + | zconsUsingSMD->setID(ZCONSUSINGSMD_ID); |
1755 | + | |
1756 | + | if (globals->haveZConsUsingSMD()){ |
1757 | + | zconsUsingSMD->setData(globals->getZconsUsingSMD()); |
1758 | + | theInfo.addProperty(zconsUsingSMD); |
1759 | + | } |
1760 | + | |
1761 | //Determine the name of ouput file and add it into SimInfo's property list | |
1762 | //Be careful, do not use inFileName, since it is a pointer which | |
1763 | //point to a string at master node, and slave nodes do not contain that string | |
# | Line 1654 | Line 1787 | void SimSetup::setupZConstraint(SimInfo& theInfo){ | |
1787 | tempParaItem.zPos = zconStamp[i]->getZpos(); | |
1788 | tempParaItem.zconsIndex = zconStamp[i]->getMolIndex(); | |
1789 | tempParaItem.kRatio = zconStamp[i]->getKratio(); | |
1790 | < | |
1790 | > | tempParaItem.havingCantVel = zconStamp[i]->haveCantVel(); |
1791 | > | tempParaItem.cantVel = zconStamp[i]->getCantVel(); |
1792 | zconsParaData->addItem(tempParaItem); | |
1793 | } | |
1794 | ||
# | Line 1672 | Line 1806 | void SimSetup::setupZConstraint(SimInfo& theInfo){ | |
1806 | //push data into siminfo, therefore, we can retrieve later | |
1807 | theInfo.addProperty(zconsParaData); | |
1808 | } | |
1809 | + | |
1810 | + | void SimSetup::makeMinimizer(){ |
1811 | + | |
1812 | + | OOPSEMinimizer* myOOPSEMinimizer; |
1813 | + | MinimizerParameterSet* param; |
1814 | + | char minimizerName[100]; |
1815 | + | |
1816 | + | for (int i = 0; i < nInfo; i++){ |
1817 | + | |
1818 | + | //prepare parameter set for minimizer |
1819 | + | param = new MinimizerParameterSet(); |
1820 | + | param->setDefaultParameter(); |
1821 | + | |
1822 | + | if (globals->haveMinimizer()){ |
1823 | + | param->setFTol(globals->getMinFTol()); |
1824 | + | } |
1825 | + | |
1826 | + | if (globals->haveMinGTol()){ |
1827 | + | param->setGTol(globals->getMinGTol()); |
1828 | + | } |
1829 | + | |
1830 | + | if (globals->haveMinMaxIter()){ |
1831 | + | param->setMaxIteration(globals->getMinMaxIter()); |
1832 | + | } |
1833 | + | |
1834 | + | if (globals->haveMinWriteFrq()){ |
1835 | + | param->setMaxIteration(globals->getMinMaxIter()); |
1836 | + | } |
1837 | + | |
1838 | + | if (globals->haveMinWriteFrq()){ |
1839 | + | param->setWriteFrq(globals->getMinWriteFrq()); |
1840 | + | } |
1841 | + | |
1842 | + | if (globals->haveMinStepSize()){ |
1843 | + | param->setStepSize(globals->getMinStepSize()); |
1844 | + | } |
1845 | + | |
1846 | + | if (globals->haveMinLSMaxIter()){ |
1847 | + | param->setLineSearchMaxIteration(globals->getMinLSMaxIter()); |
1848 | + | } |
1849 | + | |
1850 | + | if (globals->haveMinLSTol()){ |
1851 | + | param->setLineSearchTol(globals->getMinLSTol()); |
1852 | + | } |
1853 | + | |
1854 | + | strcpy(minimizerName, globals->getMinimizer()); |
1855 | + | |
1856 | + | if (!strcasecmp(minimizerName, "CG")){ |
1857 | + | myOOPSEMinimizer = new PRCGMinimizer(&(info[i]), the_ff, param); |
1858 | + | } |
1859 | + | else if (!strcasecmp(minimizerName, "SD")){ |
1860 | + | //myOOPSEMinimizer = MinimizerFactory.creatMinimizer("", &(info[i]), the_ff, param); |
1861 | + | myOOPSEMinimizer = new SDMinimizer(&(info[i]), the_ff, param); |
1862 | + | } |
1863 | + | else{ |
1864 | + | sprintf(painCave.errMsg, |
1865 | + | "SimSetup error: Unrecognized Minimizer, use Conjugate Gradient \n"); |
1866 | + | painCave.isFatal = 0; |
1867 | + | simError(); |
1868 | + | |
1869 | + | myOOPSEMinimizer = new PRCGMinimizer(&(info[i]), the_ff, param); |
1870 | + | } |
1871 | + | info[i].the_integrator = myOOPSEMinimizer; |
1872 | + | |
1873 | + | //store the minimizer into simInfo |
1874 | + | info[i].the_minimizer = myOOPSEMinimizer; |
1875 | + | info[i].has_minimizer = true; |
1876 | + | } |
1877 | + | |
1878 | + | } |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |