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 673 by chuckv, Fri Aug 8 21:22:37 2003 UTC

# 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_dr;    // The distance between each of the rho points.    
34 >    int eam_nrho;  // Number of points indexed by rho
35 >    int eam_nr;    // The number of points based on r (Both Phi(r) and Rho(r)).
36 >    double eam_rcut;  // The cutoff radius for eam.
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 >      int i;
80 >
81        // check for duplicates
82        
83        if( !strcmp( info.name, name ) ){
# Line 71 | Line 89 | namespace EAM_NS{
89          simError();
90        }
91        
92 <      if( next != NULL ) next->add(info);
92 >      if( next != NULL ) next->add(info, the_eam_rvals, the_eam_rhovals, the_eam_Frhovals);
93        else{
94          next = new LinkedAtomType();
95          strcpy(next->name, info.name);
96 <        next->mass     = info.mass;
97 <        next->epslon   = info.epslon;
98 <        next->sigma    = info.sigma;
99 <        next->ident    = info.ident;
96 >        next->mass             = info.mass;
97 >        next->lattice_constant = info.lattice_constant;
98 >        next->eam_nrho         = info.eam_nrho;
99 >        next->eam_drho         = info.eam_drho;
100 >        next->eam_nr           = info.eam_nr;
101 >        next->eam_dr           = info.eam_dr;
102 >        next->eam_rcut         = info.eam_rcut;
103 >        next->eam_ident        = info.eam_ident;
104 >        next->ident            = info.ident;
105 >
106 >        next->eam_rvals    = the_eam_rvals;
107 >        next->eam_rhovals  = the_eam_rhovals;
108 >        next->eam_Frhovals = the_eam_Frhovals;
109        }
110      }
111      
# Line 87 | Line 114 | namespace EAM_NS{
114      
115      void duplicate( atomStruct &info ){
116        strcpy(info.name, name);
117 <      info.mass     = mass;
118 <      info.epslon   = epslon;
119 <      info.sigma    = sigma;
120 <      info.ident    = ident;
121 <      info.last     = 0;
117 >      info.mass             = mass;
118 >      info.lattice_constant = lattice_constant;
119 >      info.eam_nrho         = eam_nrho;
120 >      info.eam_drho         = eam_drho;
121 >      info.eam_nr           = eam_nr;
122 >      info.eam_dr           = eam_dr;
123 >      info.eam_rcut         = eam_rcut;
124 >      info.eam_ident        = eam_ident;
125 >      info.ident            = ident;
126 >      info.last             = 0;
127      }
128  
129  
# Line 99 | Line 131 | namespace EAM_NS{
131  
132      char name[15];
133      double mass;
134 <    double epslon;
135 <    double sigma;
134 >    double lattice_constant;
135 >    int eam_nrho; // Number of points indexed by rho
136 >    double eam_drho; // The distance between each of the points indexed by rho.
137 >    int eam_nr;   // The number of points based on r (Both Phi(r) and Rho(r)).
138 >    double eam_dr;   // The distance between each of the rho points.
139 >    double eam_rcut; // The cutoff radius for eam.
140 >
141 >    double *eam_rvals;    // Z of r values
142 >    double *eam_rhovals;  // rho of r values
143 >    double *eam_Frhovals; // F of rho values
144 >    int eam_ident;        // eam identity (atomic number)
145      int ident;
146      LinkedAtomType* next;
147    };
# Line 110 | Line 151 | using namespace LJ_NS;
151  
152   }
153  
154 < using namespace LJ_NS;
154 > using namespace EAM_NS;
155  
156   //****************************************************************
157   // begins the actual forcefield stuff.  
# Line 128 | Line 169 | EAM_FF::EAM_FF(){
169    headAtomType = NULL;
170    currentAtomType = NULL;
171  
172 +  // Set eamRcut to 0.0
173 +  eamRcut = 0.0;
174 +
175    // do the funtion wrapping
176    wrapMeFF( this );
177  
# Line 138 | Line 182 | EAM_FF::EAM_FF(){
182    // Init the atomStruct mpi type
183  
184    atomStruct atomProto; // mpiPrototype
185 <  int atomBC[3] = {15,3,2};  // block counts
185 >  int atomBC[3] = {15,4,6};  // block counts
186    MPI_Aint atomDspls[3];           // displacements
187    MPI_Datatype atomMbrTypes[3];    // member mpi types
188  
189    MPI_Address(&atomProto.name, &atomDspls[0]);
190    MPI_Address(&atomProto.mass, &atomDspls[1]);
191 <  MPI_Address(&atomProto.ident, &atomDspls[2]);
191 >  MPI_Address(&atomProto.eam_nrho, &atomDspls[2]);
192    
193    atomMbrTypes[0] = MPI_CHAR;
194    atomMbrTypes[1] = MPI_DOUBLE;
# Line 201 | Line 245 | EAM_FF::EAM_FF(){
245   #ifdef IS_MPI
246    }
247    
248 <  sprintf( checkPointMsg, "LJ_FF file opened sucessfully." );
248 >  sprintf( checkPointMsg, "EAM_FF file opened sucessfully." );
249    MPIcheckPoint();
250    
251   #endif // is_mpi
# Line 223 | Line 267 | void EAM_FF::initForceField( int ljMixRule ){
267   #endif // is_mpi
268   }
269  
270 < void EAM_FF::initForceField( int ljMixRule ){
270 >
271 > void EAM_FF::calcRcut( void ){
272    
273 +  #ifdef IS_MPI
274 +  double tempEamRcut = eamRcut;
275 +  MPI_Allreduce( &tempEamRcut, &eamRcut, 1, MPI_DOUBLE, MPI_MAX,
276 +                 MPI_COMM_WORLD);
277 + #endif  //is_mpi
278 +  entry_plug->setRcut(eamRcut);
279 + }
280 +
281 +
282 + void EAM_FF::initForceField( int ljMixRule ){
283    initFortran( ljMixRule, 0 );
284   }
285  
# Line 243 | Line 298 | void EAM_FF::readParams( void ){
298   #endif // is_mpi
299   }
300  
301 +
302   void EAM_FF::readParams( void ){
303  
304    atomStruct info;
# Line 251 | Line 307 | void EAM_FF::readParams( void ){
307  
308    int i;
309    int identNum;
310 +  double *eam_rvals;    // Z of r values
311 +  double *eam_rhovals;  // rho of r values
312 +  double *eam_Frhovals; // F of rho values
313 +  char eamPotFile[1000];
314    
315  
316    bigSigma = 0.0;
# Line 262 | Line 322 | void EAM_FF::readParams( void ){
322  
323      headAtomType = new LinkedAtomType;
324      
325 <    fastForward( "AtomTypes", "initializeAtoms" );
325 >    fastForward( "AtomTypes", "eam atom readParams" );
326  
327      // we are now at the AtomTypes section.
328      
# Line 282 | Line 342 | void EAM_FF::readParams( void ){
342      
343      identNum = 1;
344      // stop reading at end of file, or at next section
345 +
346      while( readLine[0] != '#' && eof_test != NULL ){
347  
348        // toss comment lines
349        if( readLine[0] != '!' ){
350          
351          // the parser returns 0 if the line was blank
352 <        if( parseAtom( readLine, lineNum, info ) ){
352 >        if( parseAtom( readLine, lineNum, info, eamPotFile ) ){
353 >          parseEAM(info,eamPotFile, &eam_rvals,
354 >                   &eam_rhovals, &eam_Frhovals);
355            info.ident = identNum;
356 <          headAtomType->add( info );;
356 >          headAtomType->add( info, eam_rvals,
357 >                             eam_rhovals,eam_Frhovals );
358            identNum++;
359          }
360        }
361        eof_test = fgets( readLine, sizeof(readLine), frcFile );
362        lineNum++;
363      }
364 +    
365 +    
366  
367   #ifdef IS_MPI
368 +  
369      
370      // send out the linked list to all the other processes
371  
372      sprintf( checkPointMsg,
373 <             "LJ_FF atom structures read successfully." );
373 >             "EAM_FF atom structures read successfully." );
374      MPIcheckPoint();
375  
376      currentAtomType = headAtomType->next; //skip the first element who is a place holder.
# Line 314 | Line 381 | void EAM_FF::readParams( void ){
381  
382        sendFrcStruct( &info, mpiAtomStructType );
383  
384 +      // We have to now broadcast the Arrays
385 +      MPI_Bcast(currentAtomType->eam_rvals,
386 +                currentAtomType->eam_nr,
387 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
388 +      MPI_Bcast(currentAtomType->eam_rhovals,
389 +                currentAtomType->eam_nr,
390 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
391 +      MPI_Bcast(currentAtomType->eam_Frhovals,
392 +                currentAtomType->eam_nrho,
393 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
394 +
395        sprintf( checkPointMsg,
396 <               "successfully sent lJ force type: \"%s\"\n",
396 >               "successfully sent EAM force type: \"%s\"\n",
397                 info.name );
398        MPIcheckPoint();
399  
# Line 334 | Line 412 | void EAM_FF::readParams( void ){
412  
413      headAtomType = new LinkedAtomType;
414      recieveFrcStruct( &info, mpiAtomStructType );
415 <    
415 >
416      while( !info.last ){
417 +      
418 +      // allocate the arrays
419  
420 +      eam_rvals    = new double[info.eam_nr];
421 +      eam_rhovals  = new double[info.eam_nr];
422 +      eam_Frhovals = new double[info.eam_nrho];
423  
424 +      // We have to now broadcast the Arrays
425 +      MPI_Bcast(eam_rvals,
426 +                info.eam_nr,
427 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
428 +      MPI_Bcast(eam_rhovals,
429 +                info.eam_nr,
430 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
431 +      MPI_Bcast(eam_Frhovals,
432 +                info.eam_nrho,
433 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
434 +      
435  
436 <      headAtomType->add( info );
436 >      headAtomType->add( info, eam_rvals, eam_rhovals, eam_Frhovals );
437        
438        MPIcheckPoint();
439  
440        recieveFrcStruct( &info, mpiAtomStructType );
441 +
442 +
443      }
444    }
445   #endif // is_mpi
# Line 353 | Line 449 | void EAM_FF::readParams( void ){
449    int isError;
450  
451    // dummy variables
452 <  int isLJ = 1;
452 >  int isLJ = 0;
453    int isDipole = 0;
454    int isSSD = 0;
455    int isGB = 0;
456 +  int isEAM= 1;
457    double dipole = 0.0;
458 +  double eamSigma = 0.0;
459 +  double eamEpslon = 0.0;
460    
461 <  currentAtomType = headAtomType;
461 >  currentAtomType = headAtomType->next;
462    while( currentAtomType != NULL ){
463      
464      if( currentAtomType->name[0] != '\0' ){
# Line 369 | Line 468 | void EAM_FF::readParams( void ){
468                   &isSSD,
469                   &isDipole,
470                   &isGB,
471 <                 &(currentAtomType->epslon),
472 <                 &(currentAtomType->sigma),
471 >                 &isEAM,
472 >                 &eamEpslon,
473 >                 &eamSigma,
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    
538  
# Line 416 | Line 550 | void EAM_FF::initializeAtoms( int nAtoms, Atom** the_a
550      }
551      
552      the_atoms[i]->setMass( currentAtomType->mass );
419 –    the_atoms[i]->setEpslon( currentAtomType->epslon );
420 –    the_atoms[i]->setSigma( currentAtomType->sigma );
553      the_atoms[i]->setIdent( currentAtomType->ident );
554 <    the_atoms[i]->setLJ();
554 >    the_atoms[i]->setEAM();
555 >    the_atoms[i]->setEamRcut( currentAtomType->eam_rcut);
556  
557 <    if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
557 >    if (eamRcut < currentAtomType->eam_rcut) eamRcut = currentAtomType->eam_rcut;
558 >
559    }
560   }
561  
562 < void LJ_FF::initializeBonds( int nBonds, Bond** BondArray,
562 > void EAM_FF::initializeBonds( int nBonds, Bond** BondArray,
563                               bond_pair* the_bonds ){
564    
565      if( nBonds ){
# Line 514 | Line 648 | int EAM_NS::parseAtom( char *lineBuffer, int lineNum,
648  
649  
650  
651 < int EAM_NS::parseAtom( char *lineBuffer, int lineNum,  atomStruct &info ){
651 > int EAM_NS::parseAtom( char *lineBuffer, int lineNum,   atomStruct &info, char *eamPotFile ){
652  
653    char* the_token;
654    
# Line 539 | Line 673 | int EAM_NS::parseAtom( char *lineBuffer, int lineNum,
673        simError();
674      }
675          
676 <    info.epslon = atof( the_token );
677 <          
678 <    if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
679 <      sprintf( painCave.errMsg,
680 <               "Error parseing AtomTypes: line %d\n", lineNum );
676 >    strcpy( eamPotFile, the_token );
677 >    return 1;
678 >  }
679 >  else return 0;
680 > }
681 >
682 > int EAM_NS::parseEAM(atomStruct &info, char *eamPotFile,
683 >                     double **eam_rvals,
684 >                     double **eam_rhovals,
685 >                     double **eam_Frhovals){
686 >  double* myEam_rvals;
687 >  double* myEam_rhovals;
688 >  double* myEam_Frhovals;
689 >
690 >  char* ffPath_env = "FORCE_PARAM_PATH";
691 >  char* ffPath;
692 >  char* the_token;
693 >  char* eam_eof_test;
694 >  FILE *eamFile;
695 >  const int BUFFERSIZE = 3000;
696 >
697 >  char temp[200];
698 >  int linenumber;
699 >  int nReadLines;
700 >  char eam_read_buffer[BUFFERSIZE];
701 >
702 >
703 >  int i,j;
704 >
705 >  linenumber = 0;
706 >
707 >  // Open eam file
708 >  eamFile = fopen( eamPotFile, "r" );
709 >  
710 >  
711 >  if( eamFile == NULL ){
712 >    
713 >      // next see if the force path enviorment variable is set
714 >    
715 >    ffPath = getenv( ffPath_env );
716 >    if( ffPath == NULL ) {
717 >      STR_DEFINE(ffPath, FRC_PATH );
718 >    }
719 >    
720 >    
721 >    strcpy( temp, ffPath );
722 >    strcat( temp, "/" );
723 >    strcat( temp, eamPotFile );
724 >    strcpy( eamPotFile, temp );
725 >    
726 >    eamFile = fopen( eamPotFile, "r" );
727 >
728 >    
729 >    
730 >    if( eamFile == NULL ){
731 >      
732 >      sprintf( painCave.errMsg,
733 >               "Error opening the EAM force parameter file: %s\n"
734 >               "Have you tried setting the FORCE_PARAM_PATH environment "
735 >               "vairable?\n",
736 >               eamPotFile );
737        painCave.isFatal = 1;
738        simError();
739      }
550 –        
551 –    info.sigma = atof( the_token );
552 –    
553 –    return 1;
740    }
741 <  else return 0;
741 >
742 >  // First line is a comment line, read and toss it....
743 >  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
744 >  linenumber++;
745 >  if(eam_eof_test == NULL){
746 >    sprintf( painCave.errMsg,
747 >             "error in reading commment in %s\n", eamPotFile);
748 >    painCave.isFatal = 1;
749 >    simError();
750 >  }
751 >
752 >
753 >
754 >  // The Second line contains atomic number, atomic mass and a lattice constant
755 >  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
756 >  linenumber++;
757 >  if(eam_eof_test == NULL){
758 >    sprintf( painCave.errMsg,
759 >             "error in reading Identifier line in %s\n", eamPotFile);
760 >    painCave.isFatal = 1;
761 >    simError();
762 >  }
763 >
764 >
765 >
766 >    
767 >  if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
768 >    sprintf( painCave.errMsg,
769 >             "Error parseing EAM ident  line in %s\n", eamPotFile );
770 >    painCave.isFatal = 1;
771 >    simError();
772 >  }
773 >  
774 >  info.eam_ident = atoi( the_token );
775 >
776 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
777 >    sprintf( painCave.errMsg,
778 >             "Error parseing EAM mass in %s\n", eamPotFile );
779 >    painCave.isFatal = 1;
780 >    simError();
781 >  }
782 >  info.mass = atof( the_token);
783 >
784 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
785 >    sprintf( painCave.errMsg,
786 >             "Error parseing EAM Lattice Constant %s\n", eamPotFile );
787 >    painCave.isFatal = 1;
788 >    simError();
789 >  }
790 >  info.lattice_constant = atof( the_token);
791 >
792 >  // Next line is nrho, drho, nr, dr and rcut
793 >  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
794 >  if(eam_eof_test == NULL){
795 >    sprintf( painCave.errMsg,
796 >             "error in reading number of points line in %s\n", eamPotFile);
797 >    painCave.isFatal = 1;
798 >    simError();
799 >  }
800 >
801 >  if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
802 >    sprintf( painCave.errMsg,
803 >             "Error parseing EAM nrho: line in %s\n", eamPotFile );
804 >    painCave.isFatal = 1;
805 >    simError();
806 >  }
807 >  
808 >  info.eam_nrho = atoi( the_token );
809 >  
810 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
811 >    sprintf( painCave.errMsg,
812 >             "Error parseing EAM drho in %s\n", eamPotFile );
813 >    painCave.isFatal = 1;
814 >    simError();
815 >  }
816 >  info.eam_drho = atof( the_token);
817 >
818 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
819 >    sprintf( painCave.errMsg,
820 >             "Error parseing EAM # r in %s\n", eamPotFile );
821 >    painCave.isFatal = 1;
822 >    simError();
823 >  }
824 >  info.eam_nr = atoi( the_token);
825 >  
826 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
827 >    sprintf( painCave.errMsg,
828 >             "Error parseing EAM dr in %s\n", eamPotFile );
829 >    painCave.isFatal = 1;
830 >    simError();
831 >  }
832 >  info.eam_dr = atof( the_token);
833 >
834 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
835 >    sprintf( painCave.errMsg,
836 >             "Error parseing EAM rcut in %s\n", eamPotFile );
837 >    painCave.isFatal = 1;
838 >    simError();
839 >  }
840 >  info.eam_rcut = atof( the_token);
841 >
842 >
843 >
844 >
845 >
846 >  // Ok now we have to allocate point arrays and read in number of points
847 >  // Index the arrays for fortran, starting at 1
848 >  myEam_Frhovals = new double[info.eam_nrho];
849 >  myEam_rvals    = new double[info.eam_nr];
850 >  myEam_rhovals  = new double[info.eam_nr];
851 >
852 >  // Parse F of rho vals.
853 >
854 >  // Assume for now that we have a complete number of lines
855 >  nReadLines = int(info.eam_nrho/5);
856 >  
857 >
858 >
859 >  for (i=0;i<nReadLines;i++){
860 >    j = i*5;
861 >
862 >    // Read next line
863 >    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
864 >    linenumber++;
865 >    if(eam_eof_test == NULL){
866 >      sprintf( painCave.errMsg,
867 >               "error in reading EAM file %s at line %d\n",
868 >               eamPotFile,linenumber);
869 >      painCave.isFatal = 1;
870 >      simError();
871 >    }
872 >    
873 >    // Parse 5 values on each line into array
874 >    // Value 1
875 >    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
876 >      sprintf( painCave.errMsg,
877 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
878 >      painCave.isFatal = 1;
879 >      simError();
880 >    }
881 >
882 >    myEam_Frhovals[j+0] = atof( the_token );
883 >    
884 >    // Value 2
885 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
886 >      sprintf( painCave.errMsg,
887 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
888 >      painCave.isFatal = 1;
889 >      simError();
890 >    }
891 >
892 >    myEam_Frhovals[j+1] = atof( the_token );
893 >    
894 >    // Value 3
895 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
896 >      sprintf( painCave.errMsg,
897 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
898 >      painCave.isFatal = 1;
899 >      simError();
900 >    }
901 >
902 >    myEam_Frhovals[j+2] = atof( the_token );
903 >    
904 >    // Value 4
905 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
906 >      sprintf( painCave.errMsg,
907 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
908 >      painCave.isFatal = 1;
909 >      simError();
910 >    }
911 >
912 >    myEam_Frhovals[j+3] = atof( the_token );
913 >
914 >    // Value 5
915 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
916 >      sprintf( painCave.errMsg,
917 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
918 >      painCave.isFatal = 1;
919 >      simError();
920 >    }
921 >
922 >    myEam_Frhovals[j+4] = atof( the_token );
923 >    
924 >  }
925 >  // Parse Z of r vals
926 >  
927 >  // Assume for now that we have a complete number of lines
928 >  nReadLines = int(info.eam_nr/5);
929 >
930 >  for (i=0;i<nReadLines;i++){
931 >    j = i*5;
932 >
933 >    // Read next line
934 >    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
935 >    linenumber++;
936 >    if(eam_eof_test == NULL){
937 >      sprintf( painCave.errMsg,
938 >               "error in reading EAM file %s at line %d\n",
939 >               eamPotFile,linenumber);
940 >      painCave.isFatal = 1;
941 >      simError();
942 >    }
943 >    
944 >    // Parse 5 values on each line into array
945 >    // Value 1
946 >    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
947 >      sprintf( painCave.errMsg,
948 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
949 >      painCave.isFatal = 1;
950 >      simError();
951 >    }
952 >    
953 >    myEam_rvals[j+0] = atof( the_token );
954 >
955 >    // Value 2
956 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
957 >      sprintf( painCave.errMsg,
958 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
959 >      painCave.isFatal = 1;
960 >      simError();
961 >    }
962 >  
963 >    myEam_rvals[j+1] = atof( the_token );
964 >
965 >    // Value 3
966 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
967 >      sprintf( painCave.errMsg,
968 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
969 >      painCave.isFatal = 1;
970 >      simError();
971 >    }
972 >  
973 >    myEam_rvals[j+2] = atof( the_token );
974 >
975 >    // Value 4
976 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
977 >      sprintf( painCave.errMsg,
978 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
979 >      painCave.isFatal = 1;
980 >      simError();
981 >    }
982 >  
983 >    myEam_rvals[j+3] = atof( the_token );
984 >
985 >    // Value 5
986 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
987 >      sprintf( painCave.errMsg,
988 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
989 >      painCave.isFatal = 1;
990 >      simError();
991 >    }
992 >  
993 >    myEam_rvals[j+4] = atof( the_token );
994 >
995 >  }
996 >  // Parse rho of r vals
997 >
998 >  // Assume for now that we have a complete number of lines
999 >
1000 >  for (i=0;i<nReadLines;i++){
1001 >    j = i*5;
1002 >
1003 >    // Read next line
1004 >    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
1005 >    linenumber++;
1006 >    if(eam_eof_test == NULL){
1007 >      sprintf( painCave.errMsg,
1008 >               "error in reading EAM file %s at line %d\n",
1009 >               eamPotFile,linenumber);
1010 >      painCave.isFatal = 1;
1011 >      simError();
1012 >    }
1013 >  
1014 >    // Parse 5 values on each line into array
1015 >    // Value 1
1016 >    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
1017 >      sprintf( painCave.errMsg,
1018 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1019 >      painCave.isFatal = 1;
1020 >      simError();
1021 >    }
1022 >  
1023 >    myEam_rhovals[j+0] = atof( the_token );
1024 >
1025 >    // Value 2
1026 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
1027 >      sprintf( painCave.errMsg,
1028 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1029 >      painCave.isFatal = 1;
1030 >      simError();
1031 >    }
1032 >  
1033 >    myEam_rhovals[j+1] = atof( the_token );
1034 >
1035 >    // Value 3
1036 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
1037 >      sprintf( painCave.errMsg,
1038 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1039 >      painCave.isFatal = 1;
1040 >      simError();
1041 >    }
1042 >  
1043 >    myEam_rhovals[j+2] = atof( the_token );
1044 >
1045 >    // Value 4
1046 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
1047 >      sprintf( painCave.errMsg,
1048 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1049 >      painCave.isFatal = 1;
1050 >      simError();
1051 >    }
1052 >  
1053 >    myEam_rhovals[j+3] = atof( the_token );
1054 >
1055 >    // Value 5
1056 >    if ( (the_token = strtok( NULL, " \n\t,;")) == NULL){
1057 >      sprintf( painCave.errMsg,
1058 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1059 >      painCave.isFatal = 1;
1060 >      simError();
1061 >    }
1062 >  
1063 >    myEam_rhovals[j+4] = atof( the_token );
1064 >
1065 >  }
1066 >  *eam_rvals = myEam_rvals;
1067 >  *eam_rhovals = myEam_rhovals;
1068 >  *eam_Frhovals = myEam_Frhovals;
1069 >
1070 >  fclose(eamFile);
1071 >  return 0;
1072   }

Diff Legend

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