Unverified Commit 7d61fbc6 authored by Richard Berger's avatar Richard Berger
Browse files

Update pair adp

parent 5bdc3e9f

src/MANYBODY/pair_adp.cpp

+107 −108
Original line number Diff line number Diff line
@@ -29,11 +29,11 @@ 
#include "memory.h"

#include "error.h"

#include "utils.h"

#include "tokenizer.h"

#include "potential_file_reader.h"



using namespace LAMMPS_NS;



#define MAXLINE 1024



/* ---------------------------------------------------------------------- */



PairADP::PairADP(LAMMPS *lmp) : Pair(lmp)
@@ -96,7 +96,7 @@ PairADP::~PairADP() 
  if (setfl) {

    for (int i = 0; i < setfl->nelements; i++) delete [] setfl->elements[i];

    delete [] setfl->elements;

    delete [] setfl->mass;

    memory->destroy(setfl->mass);

    memory->destroy(setfl->frho);

    memory->destroy(setfl->rhor);

    memory->destroy(setfl->z2r);
@@ -453,7 +453,7 @@ void PairADP::coeff(int narg, char **arg) 
  if (setfl) {

    for (i = 0; i < setfl->nelements; i++) delete [] setfl->elements[i];

    delete [] setfl->elements;

    delete [] setfl->mass;

    memory->destroy(setfl->mass);

    memory->destroy(setfl->frho);

    memory->destroy(setfl->rhor);

    memory->destroy(setfl->z2r);
@@ -541,110 +541,129 @@ void PairADP::read_file(char *filename) 
{

  Setfl *file = setfl;



  // open potential file

  // read potential file

  if(comm->me == 0) {

    PotentialFileReader reader(lmp, filename, "ADP");



  int me = comm->me;

  FILE *fp;

  char line[MAXLINE];

    try {

      char * line = nullptr;



  if (me == 0) {

    fp = force->open_potential(filename);

    if (fp == NULL) {

      char str[128];

      snprintf(str,128,"Cannot open ADP potential file %s",filename);

      error->one(FLERR,str);

    }

  }

      reader.skip_line();

      reader.skip_line();

      reader.skip_line();



  // read and broadcast header

      // extract element names from nelements line

      line = reader.next_line(1);

      ValueTokenizer values(line);

      file->nelements = values.next_int();



      if (values.count() != file->nelements + 1)

        error->one(FLERR,"Incorrect element names in ADP potential file");



  int n;

  if (me == 0) {

    utils::sfgets(FLERR,line,MAXLINE,fp,filename,error);

    utils::sfgets(FLERR,line,MAXLINE,fp,filename,error);

    utils::sfgets(FLERR,line,MAXLINE,fp,filename,error);

    utils::sfgets(FLERR,line,MAXLINE,fp,filename,error);

    n = strlen(line) + 1;

      file->elements = new char*[file->nelements];

      for (int i = 0; i < file->nelements; i++) {

        const std::string word = values.next_string();

        const int n = word.length() + 1;

        file->elements[i] = new char[n];

        strcpy(file->elements[i], word.c_str());

      }

  MPI_Bcast(&n,1,MPI_INT,0,world);

  MPI_Bcast(line,n,MPI_CHAR,0,world);



  sscanf(line,"%d",&file->nelements);

  int nwords = utils::count_words(line);

  if (nwords != file->nelements + 1)

    error->all(FLERR,"Incorrect element names in ADP potential file");

      // 



  char **words = new char*[file->nelements+1];

  nwords = 0;

  strtok(line," \t\n\r\f");

  while ((words[nwords++] = strtok(NULL," \t\n\r\f"))) continue;

      line = reader.next_line(5);

      values = ValueTokenizer(line);

      file->nrho = values.next_int();

      file->drho = values.next_double();

      file->nr   = values.next_int();

      file->dr   = values.next_double();

      file->cut  = values.next_double();



      if ((file->nrho <= 0) || (file->nr <= 0) || (file->dr <= 0.0))

        error->one(FLERR,"Invalid EAM potential file");

      

      memory->create(file->mass, file->nelements, "pair:mass");

      memory->create(file->frho, file->nelements, file->nrho + 1, "pair:frho");

      memory->create(file->rhor, file->nelements, file->nr + 1, "pair:rhor");

      memory->create(file->z2r, file->nelements, file->nelements, file->nr + 1, "pair:z2r");

      memory->create(file->u2r, file->nelements, file->nelements, file->nr + 1, "pair:u2r");

      memory->create(file->w2r, file->nelements, file->nelements, file->nr + 1, "pair:w2r");



  file->elements = new char*[file->nelements];

      for (int i = 0; i < file->nelements; i++) {

    n = strlen(words[i]) + 1;

    file->elements[i] = new char[n];

    strcpy(file->elements[i],words[i]);

        line = reader.next_line(2);

        values = ValueTokenizer(line);

        values.next_int(); // ignore

        file->mass[i] = values.next_double();



        reader.next_dvector(file->nrho, &file->frho[i][1]);

        reader.next_dvector(file->nr, &file->rhor[i][1]);

      }

  delete [] words;



  if (me == 0) {

    utils::sfgets(FLERR,line,MAXLINE,fp,filename,error);

    sscanf(line,"%d %lg %d %lg %lg",

           &file->nrho,&file->drho,&file->nr,&file->dr,&file->cut);

      for (int i = 0; i < file->nelements; i++) {

        for (int j = 0; j <= i; j++) {

          reader.next_dvector(file->nr, &file->z2r[i][j][1]);

        }

      }



      for (int i = 0; i < file->nelements; i++) {

        for (int j = 0; j <= i; j++) {

          reader.next_dvector(file->nr, &file->u2r[i][j][1]);

        }

      }



      for (int i = 0; i < file->nelements; i++) {

        for (int j = 0; j <= i; j++) {

          reader.next_dvector(file->nr, &file->w2r[i][j][1]);

        }

      }

    } catch (TokenizerException & e) {

      error->one(FLERR, e.what());

    }

  }



  // broadcast potential information

  MPI_Bcast(&file->nelements, 1, MPI_INT, 0, world);



  MPI_Bcast(&file->nrho, 1, MPI_INT, 0, world);

  MPI_Bcast(&file->drho, 1, MPI_DOUBLE, 0, world);

  MPI_Bcast(&file->nr, 1, MPI_INT, 0, world);

  MPI_Bcast(&file->dr, 1, MPI_DOUBLE, 0, world);

  MPI_Bcast(&file->cut, 1, MPI_DOUBLE, 0, world);



  file->mass = new double[file->nelements];

  // allocate memory on other procs

  if (comm->me != 0) {

    file->elements = new char*[file->nelements];

    for (int i = 0; i < file->nelements; i++) file->elements[i] = nullptr;

    memory->create(file->mass, file->nelements, "pair:mass");

    memory->create(file->frho, file->nelements, file->nrho + 1, "pair:frho");

    memory->create(file->rhor, file->nelements, file->nr + 1, "pair:rhor");

  memory->create(file->z2r,file->nelements,file->nelements,file->nr+1,

                 "pair:z2r");

  memory->create(file->u2r,file->nelements,file->nelements,file->nr+1,

                 "pair:u2r");

  memory->create(file->w2r,file->nelements,file->nelements,file->nr+1,

                 "pair:w2r");



  int i,j,tmp;

  for (i = 0; i < file->nelements; i++) {

    if (me == 0) {

      utils::sfgets(FLERR,line,MAXLINE,fp,filename,error);

      sscanf(line,"%d %lg",&tmp,&file->mass[i]);

    memory->create(file->z2r, file->nelements, file->nelements, file->nr + 1, "pair:z2r");

    memory->create(file->u2r, file->nelements, file->nelements, file->nr + 1, "pair:u2r");

    memory->create(file->w2r, file->nelements, file->nelements, file->nr + 1, "pair:w2r");

  }



  // broadcast file->elements string array

  for (int i = 0; i < file->nelements; i++) {

    int n = strlen(file->elements[i]) + 1;

    MPI_Bcast(&n, 1, MPI_INT, 0, world);

    if (comm->me != 0) file->elements[i] = new char[n];

    MPI_Bcast(file->elements[i], n, MPI_CHAR, 0, world);

  }

    MPI_Bcast(&file->mass[i],1,MPI_DOUBLE,0,world);



    if (me == 0) grab(fp,filename,file->nrho,&file->frho[i][1]);

  // broadcast file->mass, frho, rhor

  for (int i = 0; i < file->nelements; i++) {

    MPI_Bcast(&file->mass[i],1,MPI_DOUBLE,0,world);

    MPI_Bcast(&file->frho[i][1],file->nrho,MPI_DOUBLE,0,world);

    if (me == 0) grab(fp,filename,file->nr,&file->rhor[i][1]);

    MPI_Bcast(&file->rhor[i][1],file->nr,MPI_DOUBLE,0,world);

  }



  for (i = 0; i < file->nelements; i++)

    for (j = 0; j <= i; j++) {

      if (me == 0) grab(fp,filename,file->nr,&file->z2r[i][j][1]);

  // broadcast file->z2r, u2r, w2r

  for (int i = 0; i < file->nelements; i++) {

    for (int j = 0; j <= i; j++) {

      MPI_Bcast(&file->z2r[i][j][1],file->nr,MPI_DOUBLE,0,world);

    }



  for (i = 0; i < file->nelements; i++)

    for (j = 0; j <= i; j++) {

      if (me == 0) grab(fp,filename,file->nr,&file->u2r[i][j][1]);

      MPI_Bcast(&file->u2r[i][j][1],file->nr,MPI_DOUBLE,0,world);

    }



  for (i = 0; i < file->nelements; i++)

    for (j = 0; j <= i; j++) {

      if (me == 0) grab(fp,filename,file->nr,&file->w2r[i][j][1]);

      MPI_Bcast(&file->w2r[i][j][1],file->nr,MPI_DOUBLE,0,world);

    }



  // close the potential file



  if (me == 0) fclose(fp);

  }

}



/* ----------------------------------------------------------------------
@@ -912,26 +931,6 @@ void PairADP::interpolate(int n, double delta, double *f, double **spline) 
  }

}



/* ----------------------------------------------------------------------

   grab n values from file fp and put them in list

   values can be several to a line

   only called by proc 0

------------------------------------------------------------------------- */



void PairADP::grab(FILE *fp, char *filename, int n, double *list)

{

  char *ptr;

  char line[MAXLINE];



  int i = 0;

  while (i < n) {

    utils::sfgets(FLERR,line,MAXLINE,fp,filename,error);

    ptr = strtok(line," \t\n\r\f");

    if (ptr) list[i++] = atof(ptr);

    while ((ptr = strtok(NULL," \t\n\r\f"))) list[i++] = atof(ptr);

  }

}



/* ---------------------------------------------------------------------- */



int PairADP::pack_forward_comm(int n, int *list, double *buf,

src/MANYBODY/pair_adp.h

+0 −1
Original line number Diff line number Diff line
@@ -82,7 +82,6 @@ class PairADP : public Pair { 
  void allocate();

  void array2spline();

  void interpolate(int, double, double *, double **);

  void grab(FILE *, char *, int, double *);



  void read_file(char *);

  void file2array();

CRA Git | Maintained and supported by SUSTech CRA and CCSE