ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-4/src/UseTheForce/Shapes_FF.cpp
(Generate patch)

Comparing trunk/OOPSE-4/src/UseTheForce/Shapes_FF.cpp (file contents):
Revision 1672 by gezelter, Thu Oct 28 17:20:22 2004 UTC vs.
Revision 1688 by chrisfen, Fri Oct 29 22:28:12 2004 UTC

# Line 53 | Line 53 | void Shapes_FF::calcRcut( void ){
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   }
# Line 170 | Line 172 | void Shapes_FF::readParams( void ){
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();
# Line 257 | Line 259 | void Shapes_FF::readParams( void ){
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],
# Line 265 | Line 268 | void Shapes_FF::readParams( void ){
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",
# Line 276 | Line 280 | void Shapes_FF::readParams( void ){
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" );
# Line 298 | Line 311 | void Shapes_FF::initializeAtoms( int nAtoms, Atom** th
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();
# Line 325 | Line 340 | void Shapes_FF::initializeAtoms( int nAtoms, Atom** th
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());
# Line 450 | Line 467 | void Shapes_FF::parseShapeFile(string shapeFileName, S
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));
# Line 477 | Line 494 | void Shapes_FF::parseShapeFile(string shapeFileName, S
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) {
# Line 487 | Line 503 | void Shapes_FF::parseShapeFile(string shapeFileName, S
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();
# Line 510 | Line 526 | void Shapes_FF::parseShapeFile(string shapeFileName, S
526    }
527  
528    // pass contact functions to ShapeType
529 +
530    st->setContactFuncs(functionVector);
531  
532    // now grab the range functions
# Line 532 | Line 549 | void Shapes_FF::parseShapeFile(string shapeFileName, S
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);
# Line 578 | Line 596 | void Shapes_FF::parseShapeFile(string shapeFileName, S
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);
# Line 611 | Line 630 | double Shapes_FF::findLargestContactDistance(ShapeAtom
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      
# Line 631 | Line 650 | double Shapes_FF::findLargestContactDistance(ShapeAtom
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 + }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines