53 |
|
void Shapes_FF::calcRcut( void ){ |
54 |
|
|
55 |
|
#ifdef IS_MPI |
56 |
< |
double tempShapesRcut = shapesRcut; |
56 |
> |
double tempShapesRcut = bigContact; |
57 |
|
MPI_Allreduce( &tempShapesRcut, &shapesRcut, 1, MPI_DOUBLE, MPI_MAX, |
58 |
|
MPI_COMM_WORLD); |
59 |
+ |
#else |
60 |
+ |
shapesRcut = bigContact; |
61 |
|
#endif //is_mpi |
62 |
|
entry_plug->setDefaultRcut(shapesRcut); |
63 |
|
} |
172 |
|
st = new ShapeAtomType(); |
173 |
|
|
174 |
|
st->setName(nameString); |
175 |
< |
myATID = atomTypeMap.size(); |
175 |
> |
myATID = atomTypeMap.size() + 1; |
176 |
|
st->setIdent(myATID); |
177 |
|
parseShapeFile(shapeFileName, st); |
178 |
|
st->complete(); |
259 |
|
|
260 |
|
isError = 0; |
261 |
|
myATID = at->getIdent(); |
262 |
+ |
|
263 |
|
makeShape( &nContact, &contactL[0], &contactM[0], &contactFunc[0], |
264 |
|
&contactCoeff[0], |
265 |
|
&nRange, &rangeL[0], &rangeM[0], &rangeFunc[0], |
268 |
|
&strengthFunc[0], &strengthCoeff[0], |
269 |
|
&myATID, |
270 |
|
&isError); |
271 |
+ |
|
272 |
|
if( isError ){ |
273 |
|
sprintf( painCave.errMsg, |
274 |
|
"Error initializing the \"%s\" shape in fortran\n", |
280 |
|
} |
281 |
|
} |
282 |
|
|
283 |
+ |
isError = 0; |
284 |
+ |
completeShapeFF(&isError); |
285 |
+ |
if( isError ){ |
286 |
+ |
sprintf( painCave.errMsg, |
287 |
+ |
"Error completing Shape FF in fortran\n"); |
288 |
+ |
painCave.isFatal = 1; |
289 |
+ |
simError(); |
290 |
+ |
} |
291 |
+ |
|
292 |
|
#ifdef IS_MPI |
293 |
|
sprintf( checkPointMsg, |
294 |
|
"Shapes_FF atom structures successfully sent to fortran\n" ); |
311 |
|
AtomType* at; |
312 |
|
DirectionalAtomType* dat; |
313 |
|
ShapeAtomType* sat; |
314 |
< |
double sigma; |
314 |
> |
double longCutoff; |
315 |
|
double ji[3]; |
316 |
|
double inertialMat[3][3]; |
317 |
|
Mat3x3d momInt; |
318 |
|
string myTypeString; |
319 |
|
|
320 |
+ |
bigContact = 0.0; |
321 |
+ |
|
322 |
|
for( i=0; i<nAtoms; i++ ){ |
323 |
|
|
324 |
|
myTypeString = the_atoms[i]->getType(); |
340 |
|
if ( at->isShape() ) { |
341 |
|
|
342 |
|
sat = (ShapeAtomType*)at; |
343 |
< |
sigma = findLargestContactDistance(sat); |
344 |
< |
if (sigma > bigSigma) bigSigma = sigma; |
345 |
< |
|
343 |
> |
longCutoff = findCutoffDistance(sat); |
344 |
> |
if (longCutoff > bigContact) bigContact = longCutoff; |
345 |
> |
cout << bigContact << " is the cutoff value\n"; |
346 |
> |
|
347 |
> |
entry_plug->useShapes = 1; |
348 |
|
} |
349 |
|
|
350 |
|
the_atoms[i]->setHasCharge(at->isCharge()); |
467 |
|
painCave.severity = OOPSE_ERROR; |
468 |
|
painCave.isFatal = 1; |
469 |
|
simError(); |
470 |
< |
|
470 |
> |
} else { |
471 |
|
token = strtok(inLine, delim); |
472 |
|
token = strtok(NULL, delim); |
473 |
|
st->setMass(atof(token)); |
494 |
|
if( inLine[0] != '!' && inLine[0] != '#' ){ |
495 |
|
// end marks section completion |
496 |
|
if (isEndLine(inLine)) break; |
480 |
– |
|
497 |
|
nTokens = countTokens(inLine, delim); |
498 |
|
if (nTokens != 0) { |
499 |
|
if (nTokens < 4) { |
503 |
|
painCave.severity = OOPSE_ERROR; |
504 |
|
painCave.isFatal = 1; |
505 |
|
simError(); |
506 |
< |
|
506 |
> |
} else { |
507 |
|
// read in a spherical harmonic function |
508 |
|
token = strtok(inLine, delim); |
509 |
|
rsh = new RealSphericalHarmonic(); |
526 |
|
} |
527 |
|
|
528 |
|
// pass contact functions to ShapeType |
529 |
+ |
|
530 |
|
st->setContactFuncs(functionVector); |
531 |
|
|
532 |
|
// now grab the range functions |
549 |
|
painCave.severity = OOPSE_ERROR; |
550 |
|
painCave.isFatal = 1; |
551 |
|
simError(); |
552 |
+ |
} else { |
553 |
|
|
554 |
|
// read in a spherical harmonic function |
555 |
|
token = strtok(inLine, delim); |
596 |
|
painCave.severity = OOPSE_ERROR; |
597 |
|
painCave.isFatal = 1; |
598 |
|
simError(); |
599 |
+ |
} else { |
600 |
|
|
601 |
|
// read in a spherical harmonic function |
602 |
|
token = strtok(inLine, delim); |
630 |
|
int i, j, nSteps; |
631 |
|
double theta, thetaStep, thetaMin, costheta; |
632 |
|
double phi, phiStep; |
633 |
< |
double sigma, bigSigma; |
634 |
< |
|
633 |
> |
double sigma, bs; |
634 |
> |
|
635 |
|
nSteps = 16; |
636 |
|
|
637 |
|
thetaStep = M_PI / nSteps; |
638 |
|
thetaMin = thetaStep / 2.0; |
639 |
|
phiStep = thetaStep * 2.0; |
640 |
< |
bigSigma = 0.0; |
640 |
> |
bs = 0.0; |
641 |
|
|
642 |
|
for (i = 0; i < nSteps; i++) { |
643 |
|
|
650 |
|
|
651 |
|
sigma = st->getContactValueAt(costheta, phi); |
652 |
|
|
653 |
< |
if (sigma > bigSigma) bigSigma = sigma; |
653 |
> |
if (sigma > bs) bs = sigma; |
654 |
|
} |
655 |
|
} |
656 |
|
|
657 |
< |
return bigSigma; |
657 |
> |
return bs; |
658 |
|
} |
659 |
+ |
|
660 |
+ |
|
661 |
+ |
double Shapes_FF::findCutoffDistance(ShapeAtomType* st) { |
662 |
+ |
int i, j, nSteps; |
663 |
+ |
double theta, thetaStep, thetaMin, costheta; |
664 |
+ |
double phi, phiStep; |
665 |
+ |
double sigma, range; |
666 |
+ |
double bigCut, tempCut; |
667 |
+ |
|
668 |
+ |
nSteps = 16; |
669 |
+ |
|
670 |
+ |
thetaStep = M_PI / nSteps; |
671 |
+ |
thetaMin = thetaStep / 2.0; |
672 |
+ |
phiStep = thetaStep * 2.0; |
673 |
+ |
bigCut = 0.0; |
674 |
+ |
|
675 |
+ |
for (i = 0; i < nSteps; i++) { |
676 |
+ |
|
677 |
+ |
theta = thetaMin + i * thetaStep; |
678 |
+ |
costheta = cos(theta); |
679 |
+ |
|
680 |
+ |
for (j = 0; j < nSteps; j++) { |
681 |
+ |
|
682 |
+ |
phi = j*phiStep; |
683 |
+ |
|
684 |
+ |
sigma = st->getContactValueAt(costheta, phi); |
685 |
+ |
range = st->getRangeValueAt(costheta, phi); |
686 |
+ |
|
687 |
+ |
// cutoff for a shape is taken to be (2.5*rangeVal + contactVal) |
688 |
+ |
tempCut = 2.5*range + sigma; |
689 |
+ |
|
690 |
+ |
if (tempCut > bigCut) bigCut = tempCut; |
691 |
+ |
} |
692 |
+ |
} |
693 |
+ |
|
694 |
+ |
return bigCut; |
695 |
+ |
} |