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 657 by chuckv, Wed Jul 30 21:17:01 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 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,4,6};  // 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 +
268 + void EAM_FF::calcRcut( void ){
269 +  double tempEamRcut;
270 +
271 +  entry_plug->setRcut(eamRcut);
272 + }
273 +
274 +
275   void EAM_FF::initForceField( int ljMixRule ){
227 –  
276    initFortran( ljMixRule, 0 );
277   }
278  
# Line 243 | Line 291 | void EAM_FF::readParams( void ){
291   #endif // is_mpi
292   }
293  
294 +
295   void EAM_FF::readParams( void ){
296  
297    atomStruct info;
# Line 251 | Line 300 | void EAM_FF::readParams( void ){
300  
301    int i;
302    int identNum;
303 +  double *eam_rvals;    // Z of r values
304 +  double *eam_rhovals;  // rho of r values
305 +  double *eam_Frhovals; // F of rho values
306 +  char eamPotFile[1000];
307    
308  
309    bigSigma = 0.0;
# Line 262 | Line 315 | void EAM_FF::readParams( void ){
315  
316      headAtomType = new LinkedAtomType;
317      
318 <    fastForward( "AtomTypes", "initializeAtoms" );
318 >    fastForward( "AtomTypes", "eam atom readParams" );
319  
320      // we are now at the AtomTypes section.
321      
# Line 282 | Line 335 | void EAM_FF::readParams( void ){
335      
336      identNum = 1;
337      // stop reading at end of file, or at next section
338 +
339      while( readLine[0] != '#' && eof_test != NULL ){
340  
341        // toss comment lines
342        if( readLine[0] != '!' ){
343          
344          // the parser returns 0 if the line was blank
345 <        if( parseAtom( readLine, lineNum, info ) ){
345 >        if( parseAtom( readLine, lineNum, info, eamPotFile ) ){
346 >          parseEAM(info,eamPotFile, &eam_rvals,
347 >                   &eam_rhovals, &eam_Frhovals);
348            info.ident = identNum;
349 <          headAtomType->add( info );;
349 >          headAtomType->add( info, eam_rvals,
350 >                             eam_rhovals,eam_Frhovals );
351            identNum++;
352          }
353        }
354        eof_test = fgets( readLine, sizeof(readLine), frcFile );
355        lineNum++;
356      }
357 +    
358 +    
359  
360   #ifdef IS_MPI
361 +  
362      
363      // send out the linked list to all the other processes
364  
365      sprintf( checkPointMsg,
366 <             "LJ_FF atom structures read successfully." );
366 >             "EAM_FF atom structures read successfully." );
367      MPIcheckPoint();
368  
369      currentAtomType = headAtomType->next; //skip the first element who is a place holder.
# Line 314 | Line 374 | void EAM_FF::readParams( void ){
374  
375        sendFrcStruct( &info, mpiAtomStructType );
376  
377 +      // We have to now broadcast the Arrays
378 +      MPI_Bcast(currentAtomType->eam_rvals,
379 +                currentAtomType->eam_nr,
380 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
381 +      MPI_Bcast(currentAtomType->eam_rhovals,
382 +                currentAtomType->eam_nr,
383 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
384 +      MPI_Bcast(currentAtomType->eam_Frhovals,
385 +                currentAtomType->eam_nrho,
386 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
387 +
388        sprintf( checkPointMsg,
389 <               "successfully sent lJ force type: \"%s\"\n",
389 >               "successfully sent EAM force type: \"%s\"\n",
390                 info.name );
391        MPIcheckPoint();
392  
# Line 334 | Line 405 | void EAM_FF::readParams( void ){
405  
406      headAtomType = new LinkedAtomType;
407      recieveFrcStruct( &info, mpiAtomStructType );
408 <    
408 >
409      while( !info.last ){
410 +      
411 +      // allocate the arrays
412  
413 +      eam_rvals    = new double[info.eam_nr];
414 +      eam_rhovals  = new double[info.eam_nr];
415 +      eam_Frhovals = new double[info.eam_nrho];
416  
417 +      // We have to now broadcast the Arrays
418 +      MPI_Bcast(eam_rvals,
419 +                info.eam_nr,
420 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
421 +      MPI_Bcast(eam_rhovals,
422 +                info.eam_nr,
423 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
424 +      MPI_Bcast(eam_Frhovals,
425 +                info.eam_nrho,
426 +                MPI_DOUBLE,0,MPI_COMM_WORLD);
427 +      
428  
429 <      headAtomType->add( info );
429 >      headAtomType->add( info, eam_rvals, eam_rhovals, eam_Frhovals );
430        
431        MPIcheckPoint();
432  
433        recieveFrcStruct( &info, mpiAtomStructType );
434 +
435 +
436      }
437    }
438   #endif // is_mpi
# Line 353 | Line 442 | void EAM_FF::readParams( void ){
442    int isError;
443  
444    // dummy variables
445 <  int isLJ = 1;
445 >  int isLJ = 0;
446    int isDipole = 0;
447    int isSSD = 0;
448    int isGB = 0;
449 +  int isEAM= 1;
450    double dipole = 0.0;
451 +  double eamSigma = 0.0;
452 +  double eamEpslon = 0.0;
453    
454 <  currentAtomType = headAtomType;
454 >  currentAtomType = headAtomType->next;
455    while( currentAtomType != NULL ){
456      
457      if( currentAtomType->name[0] != '\0' ){
# Line 369 | Line 461 | void EAM_FF::readParams( void ){
461                   &isSSD,
462                   &isDipole,
463                   &isGB,
464 <                 &(currentAtomType->epslon),
465 <                 &(currentAtomType->sigma),
464 >                 &isEAM,
465 >                 &eamEpslon,
466 >                 &eamSigma,
467                   &dipole,
468                   &isError );
469        if( isError ){
# Line 384 | Line 477 | void EAM_FF::readParams( void ){
477      currentAtomType = currentAtomType->next;
478    }
479        
480 <  entry_plug->useLJ = 1;
480 >  entry_plug->useLJ = 0;
481  
482 +  // Walk down again and send out EAM type
483 +  currentAtomType = headAtomType->next;
484 +  while( currentAtomType != NULL ){
485 +    
486 +    if( currentAtomType->name[0] != '\0' ){
487 +      isError = 0;
488 +      cerr << "Calling newEAMtype for type "<<currentAtomType->eam_ident <<"\n";
489 +      newEAMtype( &(currentAtomType->lattice_constant),
490 +                  &(currentAtomType->eam_nrho),
491 +                  &(currentAtomType->eam_drho),
492 +                  &(currentAtomType->eam_nr),
493 +                  &(currentAtomType->eam_dr),
494 +                  &(currentAtomType->eam_rcut),
495 +                  currentAtomType->eam_rvals,
496 +                  currentAtomType->eam_rhovals,
497 +                  currentAtomType->eam_Frhovals,
498 +                  &(currentAtomType->eam_ident),
499 +                  &isError);
500 +      cerr << "Returned from newEAMtype\n";
501 +      if( isError ){
502 +        sprintf( painCave.errMsg,
503 +                 "Error initializing the \"%s\" atom type in fortran EAM\n",
504 +                 currentAtomType->name );
505 +        painCave.isFatal = 1;
506 +        simError();
507 +      }
508 +    }
509 +    currentAtomType = currentAtomType->next;
510 +  }
511 +
512 +
513 +
514   #ifdef IS_MPI
515    sprintf( checkPointMsg,
516 <           "LJ_FF atom structures successfully sent to fortran\n" );
516 >           "EAM_FF atom structures successfully sent to fortran\n" );
517    MPIcheckPoint();
518   #endif // is_mpi
519  
520 +  cerr << "Done sending eamtypes to fortran\n";
521 +
522   }
523  
524  
# Line 416 | Line 543 | void EAM_FF::initializeAtoms( int nAtoms, Atom** the_a
543      }
544      
545      the_atoms[i]->setMass( currentAtomType->mass );
419 –    the_atoms[i]->setEpslon( currentAtomType->epslon );
420 –    the_atoms[i]->setSigma( currentAtomType->sigma );
546      the_atoms[i]->setIdent( currentAtomType->ident );
547 <    the_atoms[i]->setLJ();
547 >    the_atoms[i]->setEAM();
548  
424 –    if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
549    }
550   }
551  
552 < void LJ_FF::initializeBonds( int nBonds, Bond** BondArray,
552 > void EAM_FF::initializeBonds( int nBonds, Bond** BondArray,
553                               bond_pair* the_bonds ){
554    
555      if( nBonds ){
# Line 514 | Line 638 | int EAM_NS::parseAtom( char *lineBuffer, int lineNum,
638  
639  
640  
641 < int EAM_NS::parseAtom( char *lineBuffer, int lineNum,  atomStruct &info ){
641 > int EAM_NS::parseAtom( char *lineBuffer, int lineNum,   atomStruct &info, char *eamPotFile ){
642  
643    char* the_token;
644    
# Line 539 | Line 663 | int EAM_NS::parseAtom( char *lineBuffer, int lineNum,
663        simError();
664      }
665          
666 <    info.epslon = atof( the_token );
667 <          
668 <    if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
669 <      sprintf( painCave.errMsg,
670 <               "Error parseing AtomTypes: line %d\n", lineNum );
666 >    strcpy( eamPotFile, the_token );
667 >    return 1;
668 >  }
669 >  else return 0;
670 > }
671 >
672 > int EAM_NS::parseEAM(atomStruct &info, char *eamPotFile,
673 >                     double **eam_rvals,
674 >                     double **eam_rhovals,
675 >                     double **eam_Frhovals){
676 >  double* myEam_rvals;
677 >  double* myEam_rhovals;
678 >  double* myEam_Frhovals;
679 >
680 >  char* ffPath_env = "FORCE_PARAM_PATH";
681 >  char* ffPath;
682 >  char* the_token;
683 >  char* eam_eof_test;
684 >  FILE *eamFile;
685 >  const int BUFFERSIZE = 3000;
686 >
687 >  char temp[200];
688 >  int linenumber;
689 >  int nReadLines;
690 >  char eam_read_buffer[BUFFERSIZE];
691 >
692 >
693 >  int i,j;
694 >
695 >  linenumber = 0;
696 >
697 >  // Open eam file
698 >  eamFile = fopen( eamPotFile, "r" );
699 >  
700 >  
701 >  if( eamFile == NULL ){
702 >    
703 >      // next see if the force path enviorment variable is set
704 >    
705 >    ffPath = getenv( ffPath_env );
706 >    if( ffPath == NULL ) {
707 >      STR_DEFINE(ffPath, FRC_PATH );
708 >    }
709 >    
710 >    
711 >    strcpy( temp, ffPath );
712 >    strcat( temp, "/" );
713 >    strcat( temp, eamPotFile );
714 >    strcpy( eamPotFile, temp );
715 >    
716 >    eamFile = fopen( eamPotFile, "r" );
717 >
718 >    
719 >    
720 >    if( eamFile == NULL ){
721 >      
722 >      sprintf( painCave.errMsg,
723 >               "Error opening the EAM force parameter file: %s\n"
724 >               "Have you tried setting the FORCE_PARAM_PATH environment "
725 >               "vairable?\n",
726 >               eamPotFile );
727        painCave.isFatal = 1;
728        simError();
729      }
730 <        
731 <    info.sigma = atof( the_token );
730 >  }
731 >
732 >  // First line is a comment line, read and toss it....
733 >  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
734 >  linenumber++;
735 >  if(eam_eof_test == NULL){
736 >    sprintf( painCave.errMsg,
737 >             "error in reading commment in %s\n", eamPotFile);
738 >    painCave.isFatal = 1;
739 >    simError();
740 >  }
741 >
742 >
743 >
744 >  // The Second line contains atomic number, atomic mass and a lattice constant
745 >  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
746 >  linenumber++;
747 >  if(eam_eof_test == NULL){
748 >    sprintf( painCave.errMsg,
749 >             "error in reading Identifier line in %s\n", eamPotFile);
750 >    painCave.isFatal = 1;
751 >    simError();
752 >  }
753 >
754 >
755 >
756      
757 <    return 1;
757 >  if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
758 >    sprintf( painCave.errMsg,
759 >             "Error parseing EAM ident  line in %s\n", eamPotFile );
760 >    painCave.isFatal = 1;
761 >    simError();
762    }
763 <  else return 0;
763 >  
764 >  info.eam_ident = atoi( the_token );
765 >
766 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
767 >    sprintf( painCave.errMsg,
768 >             "Error parseing EAM mass in %s\n", eamPotFile );
769 >    painCave.isFatal = 1;
770 >    simError();
771 >  }
772 >  info.mass = atof( the_token);
773 >
774 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
775 >    sprintf( painCave.errMsg,
776 >             "Error parseing EAM Lattice Constant %s\n", eamPotFile );
777 >    painCave.isFatal = 1;
778 >    simError();
779 >  }
780 >  info.lattice_constant = atof( the_token);
781 >
782 >  // Next line is nrho, drho, nr, dr and rcut
783 >  eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
784 >  if(eam_eof_test == NULL){
785 >    sprintf( painCave.errMsg,
786 >             "error in reading number of points line in %s\n", eamPotFile);
787 >    painCave.isFatal = 1;
788 >    simError();
789 >  }
790 >
791 >  if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
792 >    sprintf( painCave.errMsg,
793 >             "Error parseing EAM nrho: line in %s\n", eamPotFile );
794 >    painCave.isFatal = 1;
795 >    simError();
796 >  }
797 >  
798 >  info.eam_nrho = atoi( the_token );
799 >  
800 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
801 >    sprintf( painCave.errMsg,
802 >             "Error parseing EAM drho in %s\n", eamPotFile );
803 >    painCave.isFatal = 1;
804 >    simError();
805 >  }
806 >  info.eam_drho = atof( the_token);
807 >
808 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
809 >    sprintf( painCave.errMsg,
810 >             "Error parseing EAM # r in %s\n", eamPotFile );
811 >    painCave.isFatal = 1;
812 >    simError();
813 >  }
814 >  info.eam_nr = atoi( the_token);
815 >  
816 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
817 >    sprintf( painCave.errMsg,
818 >             "Error parseing EAM dr in %s\n", eamPotFile );
819 >    painCave.isFatal = 1;
820 >    simError();
821 >  }
822 >  info.eam_dr = atof( the_token);
823 >
824 >  if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
825 >    sprintf( painCave.errMsg,
826 >             "Error parseing EAM rcut in %s\n", eamPotFile );
827 >    painCave.isFatal = 1;
828 >    simError();
829 >  }
830 >  info.eam_rcut = atof( the_token);
831 >
832 >
833 >
834 >
835 >
836 >  // Ok now we have to allocate point arrays and read in number of points
837 >  // Index the arrays for fortran, starting at 1
838 >  myEam_Frhovals = new double[info.eam_nrho];
839 >  myEam_rvals    = new double[info.eam_nr];
840 >  myEam_rhovals  = new double[info.eam_nr];
841 >
842 >  // Parse F of rho vals.
843 >
844 >  // Assume for now that we have a complete number of lines
845 >  nReadLines = int(info.eam_nrho/5);
846 >  
847 >
848 >
849 >  for (i=0;i<nReadLines;i++){
850 >    j = i*5;
851 >
852 >    // Read next line
853 >    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
854 >    linenumber++;
855 >    if(eam_eof_test == NULL){
856 >      sprintf( painCave.errMsg,
857 >               "error in reading EAM file %s at line %d\n",
858 >               eamPotFile,linenumber);
859 >      painCave.isFatal = 1;
860 >      simError();
861 >    }
862 >    
863 >    // Parse 5 values on each line into array
864 >    // Value 1
865 >    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
866 >      sprintf( painCave.errMsg,
867 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
868 >      painCave.isFatal = 1;
869 >      simError();
870 >    }
871 >
872 >    myEam_Frhovals[j+0] = atof( the_token );
873 >    
874 >    // Value 2
875 >    if ( (the_token = strtok( NULL, " \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+1] = atof( the_token );
883 >    
884 >    // Value 3
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+2] = atof( the_token );
893 >    
894 >    // Value 4
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+3] = atof( the_token );
903 >
904 >    // Value 5
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+4] = atof( the_token );
913 >    
914 >  }
915 >  // Parse Z of r vals
916 >  
917 >  // Assume for now that we have a complete number of lines
918 >  nReadLines = int(info.eam_nr/5);
919 >
920 >  for (i=0;i<nReadLines;i++){
921 >    j = i*5;
922 >
923 >    // Read next line
924 >    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
925 >    linenumber++;
926 >    if(eam_eof_test == NULL){
927 >      sprintf( painCave.errMsg,
928 >               "error in reading EAM file %s at line %d\n",
929 >               eamPotFile,linenumber);
930 >      painCave.isFatal = 1;
931 >      simError();
932 >    }
933 >    
934 >    // Parse 5 values on each line into array
935 >    // Value 1
936 >    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
937 >      sprintf( painCave.errMsg,
938 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
939 >      painCave.isFatal = 1;
940 >      simError();
941 >    }
942 >    
943 >    myEam_rvals[j+0] = atof( the_token );
944 >
945 >    // Value 2
946 >    if ( (the_token = strtok( NULL, " \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+1] = atof( the_token );
954 >
955 >    // Value 3
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+2] = atof( the_token );
964 >
965 >    // Value 4
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+3] = atof( the_token );
974 >
975 >    // Value 5
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+4] = atof( the_token );
984 >
985 >  }
986 >  // Parse rho of r vals
987 >
988 >  // Assume for now that we have a complete number of lines
989 >
990 >  for (i=0;i<nReadLines;i++){
991 >    j = i*5;
992 >
993 >    // Read next line
994 >    eam_eof_test = fgets(eam_read_buffer, sizeof(eam_read_buffer),eamFile);
995 >    linenumber++;
996 >    if(eam_eof_test == NULL){
997 >      sprintf( painCave.errMsg,
998 >               "error in reading EAM file %s at line %d\n",
999 >               eamPotFile,linenumber);
1000 >      painCave.isFatal = 1;
1001 >      simError();
1002 >    }
1003 >  
1004 >    // Parse 5 values on each line into array
1005 >    // Value 1
1006 >    if ( (the_token = strtok( eam_read_buffer, " \n\t,;")) == NULL){
1007 >      sprintf( painCave.errMsg,
1008 >               "Error parseing EAM nrho: line in %s\n", eamPotFile );
1009 >      painCave.isFatal = 1;
1010 >      simError();
1011 >    }
1012 >  
1013 >    myEam_rhovals[j+0] = atof( the_token );
1014 >
1015 >    // Value 2
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+1] = atof( the_token );
1024 >
1025 >    // Value 3
1026 >    if ( (the_token = strtok( eam_read_buffer, " \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+2] = atof( the_token );
1034 >
1035 >    // Value 4
1036 >    if ( (the_token = strtok( eam_read_buffer, " \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+3] = atof( the_token );
1044 >
1045 >    // Value 5
1046 >    if ( (the_token = strtok( eam_read_buffer, " \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+4] = atof( the_token );
1054 >
1055 >  }
1056 >  *eam_rvals = myEam_rvals;
1057 >  *eam_rhovals = myEam_rhovals;
1058 >  *eam_Frhovals = myEam_Frhovals;
1059 >
1060 >  fclose(eamFile);
1061 >  return 0;
1062   }

Diff Legend

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