ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/EAM_FF.cpp
(Generate patch)

Comparing trunk/OOPSE/libmdtools/EAM_FF.cpp (file contents):
Revision 499 by chuckv, Tue Apr 15 16:37:59 2003 UTC vs.
Revision 1097 by gezelter, Mon Apr 12 20:32:20 2004 UTC

# Line 1 | Line 1
1 < #include <cstdlib>
2 < #include <cstdio>
3 < #include <cstring>
1 > #include <stdlib.h>
2 > #include <stdio.h>
3 > #include <string.h>
4  
5   #include <iostream>
6   using namespace std;
# Line 28 | Line 28 | namespace EAM_NS{
28    typedef struct{
29      char name[15];
30      double mass;
31 <    double epslon;
32 <    double sigma;
31 >    double lattice_constant;
32 >    double eam_drho;  // The distance between each of the points indexed by rho.
33 >    double eam_rcut;  // The cutoff radius for eam.
34 >    double eam_dr;    // The distance between each of the rho points.    
35 >    int eam_nrho;  // Number of points indexed by rho
36 >    int eam_nr;    // The number of points based on r (Both Phi(r) and Rho(r)).
37 >    int eam_ident; // Atomic number
38      int ident;
39      int last;      //  0  -> default
40                     //  1  -> in MPI: tells nodes to stop listening
41    } atomStruct;
42  
43 <  int parseAtom( char *lineBuffer, int lineNum, atomStruct &info );
44 <  
43 >  int parseAtom( char *lineBuffer, int lineNum, atomStruct &info, char *eamPotFile );
44 >  int parseEAM( atomStruct &info, char *eamPotFile, double **eam_rvals,
45 >                double **eam_rhovals, double **eam_Frhovals);
46   #ifdef IS_MPI
47    
48    MPI_Datatype mpiAtomStructType;
# Line 48 | Line 54 | namespace EAM_NS{
54      LinkedAtomType(){
55        next = NULL;
56        name[0] = '\0';
57 +      eam_rvals    = NULL;
58 +      eam_rhovals  = NULL;
59 +      eam_Frhovals = NULL;
60      }
52 –    ~LinkedAtomType(){ if( next != NULL ) delete next; }
61  
62 +    ~LinkedAtomType(){
63 +      if( next != NULL ) delete next;
64 +      if( eam_rvals != NULL ) delete[] eam_rvals;
65 +      if( eam_rhovals != NULL ) delete[] eam_rhovals;
66 +      if( eam_Frhovals != NULL ) delete[] eam_Frhovals;
67 +    }
68 +
69      LinkedAtomType* find(char* key){
70        if( !strcmp(name, key) ) return this;
71        if( next != NULL ) return next->find(key);
# Line 58 | Line 73 | namespace EAM_NS{
73      }
74      
75  
76 <    void add( atomStruct &info ){
77 <    
76 >    void add( atomStruct &info, double *the_eam_rvals,
77 >              double *the_eam_rhovals,double *the_eam_Frhovals ){
78 >
79        // check for duplicates
80        
81        if( !strcmp( info.name, name ) ){
# Line 71 | Line 87 | namespace EAM_NS{
87          simError();
88        }
89        
90 <      if( next != NULL ) next->add(info);
90 >      if( next != NULL ) next->add(info, the_eam_rvals, the_eam_rhovals, the_eam_Frhovals);
91        else{
92          next = new LinkedAtomType();
93          strcpy(next->name, info.name);
94 <        next->mass     = info.mass;
95 <        next->epslon   = info.epslon;
96 <        next->sigma    = info.sigma;
97 <        next->ident    = info.ident;
94 >        next->mass             = info.mass;
95 >        next->lattice_constant = info.lattice_constant;
96 >        next->eam_nrho         = info.eam_nrho;
97 >        next->eam_drho         = info.eam_drho;
98 >        next->eam_nr           = info.eam_nr;
99 >        next->eam_dr           = info.eam_dr;
100 >        next->eam_rcut         = info.eam_rcut;
101 >        next->eam_ident        = info.eam_ident;
102 >        next->ident            = info.ident;
103 >
104 >        next->eam_rvals    = the_eam_rvals;
105 >        next->eam_rhovals  = the_eam_rhovals;
106 >        next->eam_Frhovals = the_eam_Frhovals;
107        }
108      }
109      
# Line 87 | Line 112 | namespace EAM_NS{
112      
113      void duplicate( atomStruct &info ){
114        strcpy(info.name, name);
115 <      info.mass     = mass;
116 <      info.epslon   = epslon;
117 <      info.sigma    = sigma;
118 <      info.ident    = ident;
119 <      info.last     = 0;
115 >      info.mass             = mass;
116 >      info.lattice_constant = lattice_constant;
117 >      info.eam_nrho         = eam_nrho;
118 >      info.eam_drho         = eam_drho;
119 >      info.eam_nr           = eam_nr;
120 >      info.eam_dr           = eam_dr;
121 >      info.eam_rcut         = eam_rcut;
122 >      info.eam_ident        = eam_ident;
123 >      info.ident            = ident;
124 >      info.last             = 0;
125      }
126  
127  
# Line 99 | Line 129 | namespace EAM_NS{
129  
130      char name[15];
131      double mass;
132 <    double epslon;
133 <    double sigma;
132 >    double lattice_constant;
133 >    int eam_nrho; // Number of points indexed by rho
134 >    double eam_drho; // The distance between each of the points indexed by rho.
135 >    int eam_nr;   // The number of points based on r (Both Phi(r) and Rho(r)).
136 >    double eam_dr;   // The distance between each of the rho points.
137 >    double eam_rcut; // The cutoff radius for eam.
138 >
139 >    double *eam_rvals;    // Z of r values
140 >    double *eam_rhovals;  // rho of r values
141 >    double *eam_Frhovals; // F of rho values
142 >    int eam_ident;        // eam identity (atomic number)
143      int ident;
144      LinkedAtomType* next;
145    };
# Line 110 | Line 149 | using namespace LJ_NS;
149  
150   }
151  
152 < using namespace LJ_NS;
152 > using namespace EAM_NS;
153  
154   //****************************************************************
155   // begins the actual forcefield stuff.  
# Line 123 | Line 162 | EAM_FF::EAM_FF(){
162    char* ffPath_env = "FORCE_PARAM_PATH";
163    char* ffPath;
164    char temp[200];
126 –  char errMsg[1000];
165  
166    headAtomType = NULL;
167    currentAtomType = NULL;
168  
169 +  // Set eamRcut to 0.0
170 +  eamRcut = 0.0;
171 +
172    // do the funtion wrapping
173    wrapMeFF( this );
174  
# Line 138 | Line 179 | EAM_FF::EAM_FF(){
179    // Init the atomStruct mpi type
180  
181    atomStruct atomProto; // mpiPrototype
182 <  int atomBC[3] = {15,3,2};  // block counts
182 >  int atomBC[3] = {15,5,5};  // block counts
183    MPI_Aint atomDspls[3];           // displacements
184    MPI_Datatype atomMbrTypes[3];    // member mpi types
185  
186    MPI_Address(&atomProto.name, &atomDspls[0]);
187    MPI_Address(&atomProto.mass, &atomDspls[1]);
188 <  MPI_Address(&atomProto.ident, &atomDspls[2]);
188 >  MPI_Address(&atomProto.eam_nrho, &atomDspls[2]);
189    
190    atomMbrTypes[0] = MPI_CHAR;
191    atomMbrTypes[1] = MPI_DOUBLE;
# Line 201 | Line 242 | EAM_FF::EAM_FF(){
242   #ifdef IS_MPI
243    }
244    
245 <  sprintf( checkPointMsg, "LJ_FF file opened sucessfully." );
245 >  sprintf( checkPointMsg, "EAM_FF file opened sucessfully." );
246    MPIcheckPoint();
247    
248   #endif // is_mpi
# Line 223 | Line 264 | void EAM_FF::initForceField( int ljMixRule ){
264   #endif // is_mpi
265   }
266  
267 < void EAM_FF::initForceField( int ljMixRule ){
267 >
268 > void EAM_FF::calcRcut( void ){
269    
270 + #ifdef IS_MPI
271 +  double tempEamRcut = eamRcut;
272 +  MPI_Allreduce( &tempEamRcut, &eamRcut, 1, MPI_DOUBLE, MPI_MAX,
273 +                 MPI_COMM_WORLD);
274 + #endif  //is_mpi
275 +  entry_plug->setDefaultRcut(eamRcut);
276 + }
277 +
278 +
279 + void EAM_FF::initForceField( int ljMixRule ){
280    initFortran( ljMixRule, 0 );
281   }
282  
# Line 243 | Line 295 | void EAM_FF::readParams( void ){
295   #endif // is_mpi
296   }
297  
298 +
299   void EAM_FF::readParams( void ){
300  
301    atomStruct info;
302    info.last = 1; // initialize last to have the last set.
303                   // if things go well, last will be set to 0
304  
252 –  int i;
305    int identNum;
306 +  double *eam_rvals;    // Z of r values
307 +  double *eam_rhovals;  // rho of r values
308 +  double *eam_Frhovals; // F of rho values
309 +  char eamPotFile[1000];
310    
311  
312    bigSigma = 0.0;
# Line 262 | Line 318 | void EAM_FF::readParams( void ){
318  
319      headAtomType = new LinkedAtomType;
320      
321 <    fastForward( "AtomTypes", "initializeAtoms" );
321 >    fastForward( "AtomTypes", "eam atom readParams" );
322  
323      // we are now at the AtomTypes section.
324      
# Line 282 | Line 338 | void EAM_FF::readParams( void ){
338      
339      identNum = 1;
340      // stop reading at end of file, or at next section
341 +
342      while( readLine[0] != '#' && eof_test != NULL ){
343  
344        // toss comment lines
345        if( readLine[0] != '!' ){
346          
347          // the parser returns 0 if the line was blank
348 <        if( parseAtom( readLine, lineNum, info ) ){
348 >        if( parseAtom( readLine, lineNum, info, eamPotFile ) ){
349 >          parseEAM(info,eamPotFile, &eam_rvals,
350 >                   &eam_rhovals, &eam_Frhovals);
351            info.ident = identNum;
352 <          headAtomType->add( info );;
352 >          headAtomType->add( info, eam_rvals,
353 >                             eam_rhovals,eam_Frhovals );
354            identNum++;
355          }
356        }
357        eof_test = fgets( readLine, sizeof(readLine), frcFile );
358        lineNum++;
359      }
360 +    
361 +    
362  
363   #ifdef IS_MPI
364 +  
365      
366      // send out the linked list to all the other processes
367  
368      sprintf( checkPointMsg,
369 <             "LJ_FF atom structures read successfully." );
369 >             "EAM_FF atom structures read successfully." );
370      MPIcheckPoint();
371  
372      currentAtomType = headAtomType->next; //skip the first element who is a place holder.
# Line 314 | Line 377 | void EAM_FF::readParams( void ){
377  
378        sendFrcStruct( &info, mpiAtomStructType );
379  
380 +      // We have to now broadcast the Arrays
381 +      MPI_Bcast(currentAtomType->eam_rvals,
382 +                currentAtomType->eam_nr,
383 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
384 +      MPI_Bcast(currentAtomType->eam_rhovals,
385 +                currentAtomType->eam_nr,
386 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
387 +      MPI_Bcast(currentAtomType->eam_Frhovals,
388 +                currentAtomType->eam_nrho,
389 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
390 +
391        sprintf( checkPointMsg,
392 <               "successfully sent lJ force type: \"%s\"\n",
392 >               "successfully sent EAM force type: \"%s\"\n",
393                 info.name );
394        MPIcheckPoint();
395  
# Line 333 | Line 407 | void EAM_FF::readParams( void ){
407      MPIcheckPoint();
408  
409      headAtomType = new LinkedAtomType;
410 <    recieveFrcStruct( &info, mpiAtomStructType );
411 <    
410 >    receiveFrcStruct( &info, mpiAtomStructType );
411 >
412      while( !info.last ){
413 +      
414 +      // allocate the arrays
415  
416 +      eam_rvals    = new double[info.eam_nr];
417 +      eam_rhovals  = new double[info.eam_nr];
418 +      eam_Frhovals = new double[info.eam_nrho];
419  
420 +      // We have to now broadcast the Arrays
421 +      MPI_Bcast(eam_rvals,
422 +                info.eam_nr,
423 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
424 +      MPI_Bcast(eam_rhovals,
425 +                info.eam_nr,
426 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
427 +      MPI_Bcast(eam_Frhovals,
428 +                info.eam_nrho,
429 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
430 +      
431  
432 <      headAtomType->add( info );
432 >      headAtomType->add( info, eam_rvals, eam_rhovals, eam_Frhovals );
433        
434        MPIcheckPoint();
435  
436 <      recieveFrcStruct( &info, mpiAtomStructType );
436 >      receiveFrcStruct( &info, mpiAtomStructType );
437 >
438 >
439      }
440    }
441   #endif // is_mpi
# Line 353 | Line 445 | void EAM_FF::readParams( void ){
445    int isError;
446  
447    // dummy variables
448 <  int isLJ = 1;
448 >  int isLJ = 0;
449    int isDipole = 0;
450    int isSSD = 0;
451    int isGB = 0;
452 +  int isEAM = 1;
453 +  int isCharge = 0;
454    double dipole = 0.0;
455 +  double charge = 0.0;
456 +  double eamSigma = 0.0;
457 +  double eamEpslon = 0.0;
458    
459 <  currentAtomType = headAtomType;
459 >  currentAtomType = headAtomType->next;
460    while( currentAtomType != NULL ){
461      
462      if( currentAtomType->name[0] != '\0' ){
# Line 369 | Line 466 | void EAM_FF::readParams( void ){
466                   &isSSD,
467                   &isDipole,
468                   &isGB,
469 <                 &(currentAtomType->epslon),
470 <                 &(currentAtomType->sigma),
469 >                 &isEAM,
470 >                 &isCharge,
471 >                 &eamEpslon,
472 >                 &eamSigma,
473 >                 &charge,
474                   &dipole,
475                   &isError );
476        if( isError ){
# Line 384 | Line 484 | void EAM_FF::readParams( void ){
484      currentAtomType = currentAtomType->next;
485    }
486        
487 <  entry_plug->useLJ = 1;
487 >  entry_plug->useLJ = 0;
488 >  entry_plug->useEAM = 1;
489 >  // Walk down again and send out EAM type
490 >  currentAtomType = headAtomType->next;
491 >  while( currentAtomType != NULL ){
492 >    
493 >    if( currentAtomType->name[0] != '\0' ){
494 >      isError = 0;
495  
496 +      newEAMtype( &(currentAtomType->lattice_constant),
497 +                  &(currentAtomType->eam_nrho),
498 +                  &(currentAtomType->eam_drho),
499 +                  &(currentAtomType->eam_nr),
500 +                  &(currentAtomType->eam_dr),
501 +                  &(currentAtomType->eam_rcut),
502 +                  currentAtomType->eam_rvals,
503 +                  currentAtomType->eam_rhovals,
504 +                  currentAtomType->eam_Frhovals,
505 +                  &(currentAtomType->eam_ident),
506 +                  &isError);
507 +
508 +      if( isError ){
509 +        sprintf( painCave.errMsg,
510 +                 "Error initializing the \"%s\" atom type in fortran EAM\n",
511 +                 currentAtomType->name );
512 +        painCave.isFatal = 1;
513 +        simError();
514 +      }
515 +    }
516 +    currentAtomType = currentAtomType->next;
517 +  }
518 +
519 +
520 +
521   #ifdef IS_MPI
522    sprintf( checkPointMsg,
523 <           "LJ_FF atom structures successfully sent to fortran\n" );
523 >           "EAM_FF atom structures successfully sent to fortran\n" );
524    MPIcheckPoint();
525   #endif // is_mpi
526  
527 +
528 +
529   }
530  
531  
532   void EAM_FF::initializeAtoms( int nAtoms, Atom** the_atoms ){
533    
534    int i;
535 <
535 >  
536    // initialize the atoms
537    
404 –
405 –  Atom* thisAtom;
406 –
538    for( i=0; i<nAtoms; i++ ){
539      
540      currentAtomType = headAtomType->find( the_atoms[i]->getType() );
# Line 416 | Line 547 | void EAM_FF::initializeAtoms( int nAtoms, Atom** the_a
547      }
548      
549      the_atoms[i]->setMass( currentAtomType->mass );
419 –    the_atoms[i]->setEpslon( currentAtomType->epslon );
420 –    the_atoms[i]->setSigma( currentAtomType->sigma );
550      the_atoms[i]->setIdent( currentAtomType->ident );
422 –    the_atoms[i]->setLJ();
551  
552 <    if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
552 >    if (eamRcut < currentAtomType->eam_rcut) eamRcut = currentAtomType->eam_rcut;
553 >    
554    }
555   }
556  
557 < void LJ_FF::initializeBonds( int nBonds, Bond** BondArray,
557 > void EAM_FF::initializeBonds( int nBonds, Bond** BondArray,
558                               bond_pair* the_bonds ){
559    
560      if( nBonds ){
# Line 514 | Line 643 | int EAM_NS::parseAtom( char *lineBuffer, int lineNum,
643  
644  
645  
646 < int EAM_NS::parseAtom( char *lineBuffer, int lineNum,  atomStruct &info ){
646 > int EAM_NS::parseAtom( char *lineBuffer, int lineNum,   atomStruct &info, char *eamPotFile ){
647  
648    char* the_token;
649    
# Line 532 | Line 661 | int EAM_NS::parseAtom( char *lineBuffer, int lineNum,
661      
662      info.mass = atof( the_token );
663      
535 –    if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
536 –      sprintf( painCave.errMsg,
537 –               "Error parseing AtomTypes: line %d\n", lineNum );
538 –      painCave.isFatal = 1;
539 –      simError();
540 –    }
541 –        
542 –    info.epslon = atof( the_token );
543 –          
664      if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
665        sprintf( painCave.errMsg,
666                 "Error parseing AtomTypes: line %d\n", lineNum );
# Line 548 | Line 668 | int EAM_NS::parseAtom( char *lineBuffer, int lineNum,
668        simError();
669      }
670          
671 <    info.sigma = atof( the_token );
552 <    
671 >    strcpy( eamPotFile, the_token );
672      return 1;
673    }
674    else return 0;
675 + }
676 +
677 + int EAM_NS::parseEAM(atomStruct &info, char *eamPotFile,
678 +                     double **eam_rvals,
679 +                     double **eam_rhovals,
680 +                     double **eam_Frhovals){
681 +  double* myEam_rvals;
682 +  double* myEam_rhovals;
683 +  double* myEam_Frhovals;
684 +
685 +  char* ffPath_env = "FORCE_PARAM_PATH";
686 +  char* ffPath;
687 +  char* the_token;
688 +  char* eam_eof_test;
689 +  FILE *eamFile;
690 +  const int BUFFERSIZE = 3000;
691 +
692 +  char temp[200];
693 +  int linenumber;
694 +  int nReadLines;
695 +  char eam_read_buffer[BUFFERSIZE];
696 +
697 +
698 +  int i,j;
699 +
700 +  linenumber = 0;
701 +
702 +  // Open eam file
703 +  eamFile = fopen( eamPotFile, "r" );
704 +  
705 +  
706 +  if( eamFile == NULL ){
707 +    
708 +      // next see if the force path enviorment variable is set
709 +    
710 +    ffPath = getenv( ffPath_env );
711 +    if( ffPath == NULL ) {
712 +      STR_DEFINE(ffPath, FRC_PATH );
713 +    }
714 +    
715 +    
716 +    strcpy( temp, ffPath );
717 +    strcat( temp, "/" );
718 +    strcat( temp, eamPotFile );
719 +    strcpy( eamPotFile, temp );
720 +    
721 +    eamFile = fopen( eamPotFile, "r" );
722 +
723 +    
724 +    
725 +    if( eamFile == NULL ){
726 +      
727 +      sprintf( painCave.errMsg,
728 +               "Error opening the EAM force parameter file: %s\n"
729 +               "Have you tried setting the FORCE_PARAM_PATH environment "
730 +               "vairable?\n",
731 +               eamPotFile );
732 +      painCave.isFatal = 1;
733 +      simError();
734 +    }
735 +  }
736 +
737 +  // First line is a comment line, read and toss it....
738 +  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
739 +  linenumber++;
740 +  if(eam_eof_test == NULL){
741 +    sprintf( painCave.errMsg,
742 +             "error in reading commment in %s\n", eamPotFile);
743 +    painCave.isFatal = 1;
744 +    simError();
745 +  }
746 +
747 +
748 +
749 +  // The Second line contains atomic number, atomic mass and a lattice constant
750 +  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
751 +  linenumber++;
752 +  if(eam_eof_test == NULL){
753 +    sprintf( painCave.errMsg,
754 +             "error in reading Identifier line in %s\n", eamPotFile);
755 +    painCave.isFatal = 1;
756 +    simError();
757 +  }
758 +
759 +
760 +
761 +    
762 +  if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
763 +    sprintf( painCave.errMsg,
764 +             "Error parseing EAM ident  line in %s\n", eamPotFile );
765 +    painCave.isFatal = 1;
766 +    simError();
767 +  }
768 +  
769 +  info.eam_ident = atoi( the_token );
770 +
771 +  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
772 +    sprintf( painCave.errMsg,
773 +             "Error parseing EAM mass in %s\n", eamPotFile );
774 +    painCave.isFatal = 1;
775 +    simError();
776 +  }
777 +  info.mass = atof( the_token);
778 +
779 +  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
780 +    sprintf( painCave.errMsg,
781 +             "Error parseing EAM Lattice Constant %s\n", eamPotFile );
782 +    painCave.isFatal = 1;
783 +    simError();
784 +  }
785 +  info.lattice_constant = atof( the_token);
786 +
787 +  // Next line is nrho, drho, nr, dr and rcut
788 +  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
789 +  if(eam_eof_test == NULL){
790 +    sprintf( painCave.errMsg,
791 +             "error in reading number of points line in %s\n", eamPotFile);
792 +    painCave.isFatal = 1;
793 +    simError();
794 +  }
795 +
796 +  if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
797 +    sprintf( painCave.errMsg,
798 +             "Error parseing EAM nrho: line in %s\n", eamPotFile );
799 +    painCave.isFatal = 1;
800 +    simError();
801 +  }
802 +  
803 +  info.eam_nrho = atoi( the_token );
804 +  
805 +  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
806 +    sprintf( painCave.errMsg,
807 +             "Error parseing EAM drho in %s\n", eamPotFile );
808 +    painCave.isFatal = 1;
809 +    simError();
810 +  }
811 +  info.eam_drho = atof( the_token);
812 +
813 +  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
814 +    sprintf( painCave.errMsg,
815 +             "Error parseing EAM # r in %s\n", eamPotFile );
816 +    painCave.isFatal = 1;
817 +    simError();
818 +  }
819 +  info.eam_nr = atoi( the_token);
820 +  
821 +  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
822 +    sprintf( painCave.errMsg,
823 +             "Error parseing EAM dr in %s\n", eamPotFile );
824 +    painCave.isFatal = 1;
825 +    simError();
826 +  }
827 +  info.eam_dr = atof( the_token);
828 +
829 +  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
830 +    sprintf( painCave.errMsg,
831 +             "Error parseing EAM rcut in %s\n", eamPotFile );
832 +    painCave.isFatal = 1;
833 +    simError();
834 +  }
835 +  info.eam_rcut = atof( the_token);
836 +
837 +
838 +
839 +
840 +
841 +  // Ok now we have to allocate point arrays and read in number of points
842 +  // Index the arrays for fortran, starting at 1
843 +  myEam_Frhovals = new double[info.eam_nrho];
844 +  myEam_rvals    = new double[info.eam_nr];
845 +  myEam_rhovals  = new double[info.eam_nr];
846 +
847 +  // Parse F of rho vals.
848 +
849 +  // Assume for now that we have a complete number of lines
850 +  nReadLines = int(info.eam_nrho/5);
851 +  
852 +
853 +
854 +  for (i=0;i<nReadLines;i++){
855 +    j = i*5;
856 +
857 +    // Read next line
858 +    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
859 +    linenumber++;
860 +    if(eam_eof_test == NULL){
861 +      sprintf( painCave.errMsg,
862 +               "error in reading EAM file %s at line %d\n",
863 +               eamPotFile,linenumber);
864 +      painCave.isFatal = 1;
865 +      simError();
866 +    }
867 +    
868 +    // Parse 5 values on each line into array
869 +    // Value 1
870 +    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
871 +      sprintf( painCave.errMsg,
872 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
873 +      painCave.isFatal = 1;
874 +      simError();
875 +    }
876 +
877 +    myEam_Frhovals[j+0] = atof( the_token );
878 +    
879 +    // Value 2
880 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
881 +      sprintf( painCave.errMsg,
882 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
883 +      painCave.isFatal = 1;
884 +      simError();
885 +    }
886 +
887 +    myEam_Frhovals[j+1] = atof( the_token );
888 +    
889 +    // Value 3
890 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
891 +      sprintf( painCave.errMsg,
892 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
893 +      painCave.isFatal = 1;
894 +      simError();
895 +    }
896 +
897 +    myEam_Frhovals[j+2] = atof( the_token );
898 +    
899 +    // Value 4
900 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
901 +      sprintf( painCave.errMsg,
902 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
903 +      painCave.isFatal = 1;
904 +      simError();
905 +    }
906 +
907 +    myEam_Frhovals[j+3] = atof( the_token );
908 +
909 +    // Value 5
910 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
911 +      sprintf( painCave.errMsg,
912 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
913 +      painCave.isFatal = 1;
914 +      simError();
915 +    }
916 +
917 +    myEam_Frhovals[j+4] = atof( the_token );
918 +    
919 +  }
920 +  // Parse Z of r vals
921 +  
922 +  // Assume for now that we have a complete number of lines
923 +  nReadLines = int(info.eam_nr/5);
924 +
925 +  for (i=0;i<nReadLines;i++){
926 +    j = i*5;
927 +
928 +    // Read next line
929 +    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
930 +    linenumber++;
931 +    if(eam_eof_test == NULL){
932 +      sprintf( painCave.errMsg,
933 +               "error in reading EAM file %s at line %d\n",
934 +               eamPotFile,linenumber);
935 +      painCave.isFatal = 1;
936 +      simError();
937 +    }
938 +    
939 +    // Parse 5 values on each line into array
940 +    // Value 1
941 +    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
942 +      sprintf( painCave.errMsg,
943 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
944 +      painCave.isFatal = 1;
945 +      simError();
946 +    }
947 +    
948 +    myEam_rvals[j+0] = atof( the_token );
949 +
950 +    // Value 2
951 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
952 +      sprintf( painCave.errMsg,
953 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
954 +      painCave.isFatal = 1;
955 +      simError();
956 +    }
957 +  
958 +    myEam_rvals[j+1] = atof( the_token );
959 +
960 +    // Value 3
961 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
962 +      sprintf( painCave.errMsg,
963 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
964 +      painCave.isFatal = 1;
965 +      simError();
966 +    }
967 +  
968 +    myEam_rvals[j+2] = atof( the_token );
969 +
970 +    // Value 4
971 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
972 +      sprintf( painCave.errMsg,
973 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
974 +      painCave.isFatal = 1;
975 +      simError();
976 +    }
977 +  
978 +    myEam_rvals[j+3] = atof( the_token );
979 +
980 +    // Value 5
981 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
982 +      sprintf( painCave.errMsg,
983 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
984 +      painCave.isFatal = 1;
985 +      simError();
986 +    }
987 +  
988 +    myEam_rvals[j+4] = atof( the_token );
989 +
990 +  }
991 +  // Parse rho of r vals
992 +
993 +  // Assume for now that we have a complete number of lines
994 +
995 +  for (i=0;i<nReadLines;i++){
996 +    j = i*5;
997 +
998 +    // Read next line
999 +    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
1000 +    linenumber++;
1001 +    if(eam_eof_test == NULL){
1002 +      sprintf( painCave.errMsg,
1003 +               "error in reading EAM file %s at line %d\n",
1004 +               eamPotFile,linenumber);
1005 +      painCave.isFatal = 1;
1006 +      simError();
1007 +    }
1008 +  
1009 +    // Parse 5 values on each line into array
1010 +    // Value 1
1011 +    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
1012 +      sprintf( painCave.errMsg,
1013 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1014 +      painCave.isFatal = 1;
1015 +      simError();
1016 +    }
1017 +  
1018 +    myEam_rhovals[j+0] = atof( the_token );
1019 +
1020 +    // Value 2
1021 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
1022 +      sprintf( painCave.errMsg,
1023 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1024 +      painCave.isFatal = 1;
1025 +      simError();
1026 +    }
1027 +  
1028 +    myEam_rhovals[j+1] = atof( the_token );
1029 +
1030 +    // Value 3
1031 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
1032 +      sprintf( painCave.errMsg,
1033 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1034 +      painCave.isFatal = 1;
1035 +      simError();
1036 +    }
1037 +  
1038 +    myEam_rhovals[j+2] = atof( the_token );
1039 +
1040 +    // Value 4
1041 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
1042 +      sprintf( painCave.errMsg,
1043 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1044 +      painCave.isFatal = 1;
1045 +      simError();
1046 +    }
1047 +  
1048 +    myEam_rhovals[j+3] = atof( the_token );
1049 +
1050 +    // Value 5
1051 +    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
1052 +      sprintf( painCave.errMsg,
1053 +               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1054 +      painCave.isFatal = 1;
1055 +      simError();
1056 +    }
1057 +  
1058 +    myEam_rhovals[j+4] = atof( the_token );
1059 +
1060 +  }
1061 +  *eam_rvals = myEam_rvals;
1062 +  *eam_rhovals = myEam_rhovals;
1063 +  *eam_Frhovals = myEam_Frhovals;
1064 +
1065 +  fclose(eamFile);
1066 +  return 0;
1067   }

Diff Legend

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