From 1ed74de9a87f4303b50cbf55604fc38a5673ea53 Mon Sep 17 00:00:00 2001 From: MaziarHeidari Date: Fri, 24 Jun 2016 16:05:27 +0200 Subject: [PATCH] HADRESS_Version_24_Jun_2016 --- src/atom.cpp | 26 ++++++++++++++++++++++++++ src/atom.h | 12 ++++++++++++ src/molecule.cpp | 45 ++++++++++++++++++++++++++++++++++++++++++++- src/molecule.h | 9 +++++++++ src/read_data.cpp | 8 +++++++- 5 files changed, 98 insertions(+), 2 deletions(-) diff --git a/src/atom.cpp b/src/atom.cpp index 8df360221c..4da03147b5 100644 --- a/src/atom.cpp +++ b/src/atom.cpp @@ -93,6 +93,15 @@ Atom::Atom(LAMMPS *lmp) : Pointers(lmp) rho = drho = e = de = cv = NULL; vest = NULL; + // USER-HAdResS package + + lambdaH = NULL; + gradlambdaH = NULL; + replambdaH = NULL; + moltypeH = NULL; + comH = NULL; + nmoltypesH = 0; + // USER-DPD uCond = uMech = uChem = uCG = uCGnew = NULL; @@ -240,6 +249,12 @@ Atom::~Atom() memory->destroy(molindex); memory->destroy(molatom); + memory->destroy(lambdaH); + memory->destroy(replambdaH); + memory->destroy(moltypeH); + memory->destroy(gradlambdaH); + memory->destroy(comH); + memory->destroy(q); memory->destroy(mu); memory->destroy(omega); @@ -392,6 +407,12 @@ void Atom::create_avec(const char *style, int narg, char **arg, int trysuffix) molecule_flag = 0; q_flag = mu_flag = 0; + + // USER-HAdResS + + replambdaH_flag = 0; + moltypeH_flag = 0; + omega_flag = torque_flag = angmom_flag = 0; radius_flag = rmass_flag = 0; ellipsoid_flag = line_flag = tri_flag = body_flag = 0; @@ -1608,6 +1629,8 @@ void Atom::add_molecule_atom(Molecule *onemol, int iatom, if (onemol->qflag && q_flag) q[ilocal] = onemol->q[iatom]; if (onemol->radiusflag && radius_flag) radius[ilocal] = onemol->radius[iatom]; if (onemol->rmassflag && rmass_flag) rmass[ilocal] = onemol->rmass[iatom]; + if (onemol->replambdaHflag && replambdaH_flag) replambdaH[ilocal] = onemol->replambdaH[iatom]; + if (onemol->moltypeHflag && moltypeH_flag) moltypeH[ilocal] = onemol->moltypeH[iatom]; else if (rmass_flag) rmass[ilocal] = 4.0*MY_PI/3.0 * radius[ilocal]*radius[ilocal]*radius[ilocal]; @@ -2082,6 +2105,9 @@ void *Atom::extract(char *name) if (strcmp(name,"f") == 0) return (void *) f; if (strcmp(name,"molecule") == 0) return (void *) molecule; if (strcmp(name,"q") == 0) return (void *) q; + if (strcmp(name,"lambdaH") == 0) return (void *) lambdaH; + if (strcmp(name,"gradlambdaH") == 0) return (void *) gradlambdaH; + if (strcmp(name,"replambdaH") == 0) return (void *) replambdaH; if (strcmp(name,"mu") == 0) return (void *) mu; if (strcmp(name,"omega") == 0) return (void *) omega; if (strcmp(name,"angmom") == 0) return (void *) angmom; diff --git a/src/atom.h b/src/atom.h index f2415e1db5..55ff0518e5 100644 --- a/src/atom.h +++ b/src/atom.h @@ -59,6 +59,14 @@ class Atom : protected Pointers { double *radius,*rmass; int *ellipsoid,*line,*tri,*body; + //USER-HAdResS package + + double *lambdaH,**gradlambdaH; + int *replambdaH; + int *moltypeH; + int nmoltypesH; + double **comH; + // PERI package double *vfrac,*s0; @@ -141,6 +149,10 @@ class Atom : protected Pointers { int rho_flag,e_flag,cv_flag,vest_flag; int dpd_flag; + // USER-HAdResS + + int replambdaH_flag, moltypeH_flag; + // USER-SMD package int smd_flag; diff --git a/src/molecule.cpp b/src/molecule.cpp index 64ebf2f4ae..8b954a97b3 100644 --- a/src/molecule.cpp +++ b/src/molecule.cpp @@ -503,7 +503,20 @@ void Molecule::read(int flag) if (flag) masses(line); else skip_lines(natoms,line); - } else if (strcmp(keyword,"Bonds") == 0) { + } else if (strcmp(keyword,"replambdaH") == 0) { + // USER-HAdResS Package + + replambdaHflag = 1; + if (flag) representative_atom(line); + else skip_lines(natoms,line); + } else if (strcmp(keyword,"moltypeH") == 0) { + // USER-HAdResS Package + + moltypeHflag = 1; + if (flag) moltype_atom(line); + else skip_lines(natoms,line); + + } else if (strcmp(keyword,"Bonds") == 0) { if (nbonds == 0) error->all(FLERR,"Molecule file has bonds but no nbonds setting"); bondflag = tag_require = 1; @@ -685,6 +698,30 @@ void Molecule::charges(char *line) } } +// USER-HAdResS Package + +void Molecule::representative_atom(char *line) +{ + int tmp; + for (int i = 0; i < natoms; i++) { + readline(line); + sscanf(line,"%d %d",&tmp,&replambdaH[i]); + printf("i=%d rep = %d",i,replambdaH[i]); + } +} + +// USER-HAdResS Package + +void Molecule::moltype_atom(char *line) +{ + int tmp; + for (int i = 0; i < natoms; i++) { + readline(line); + sscanf(line,"%d %d",&tmp,&moltypeH[i]); + printf("i=%d rep = %d",i,moltypeH[i]); + } +} + /* ---------------------------------------------------------------------- read diameters from file and set radii ------------------------------------------------------------------------- */ @@ -1342,6 +1379,8 @@ void Molecule::check_attributes(int flag) // warn if not a match int mismatch = 0; + if (onemol->replambdaHflag && !atom->replambdaH_flag) mismatch = 1; + if (onemol->moltypeHflag && !atom->moltypeH_flag) mismatch = 1; if (onemol->qflag && !atom->q_flag) mismatch = 1; if (onemol->radiusflag && !atom->radius_flag) mismatch = 1; if (onemol->rmassflag && !atom->rmass_flag) mismatch = 1; @@ -1460,6 +1499,8 @@ void Molecule::allocate() if (qflag) memory->create(q,natoms,"molecule:q"); if (radiusflag) memory->create(radius,natoms,"molecule:radius"); if (rmassflag) memory->create(rmass,natoms,"molecule:rmass"); + if (replambdaHflag) memory->create(replambdaH,natoms,"molecule:replambdaH"); + if (moltypeHflag) memory->create(moltypeH,natoms,"molecule:moltypeH"); // always allocate num_bond,num_angle,etc and special+nspecial // even if not in molecule file, initialize to 0 @@ -1546,6 +1587,8 @@ void Molecule::deallocate() memory->destroy(type); memory->destroy(q); memory->destroy(radius); + memory->destroy(replambdaH); + memory->destroy(moltypeH); memory->destroy(rmass); memory->destroy(num_bond); diff --git a/src/molecule.h b/src/molecule.h index 0bbe684636..fb7f86622e 100644 --- a/src/molecule.h +++ b/src/molecule.h @@ -42,6 +42,7 @@ class Molecule : protected Pointers { // 1 if attribute defined in file, 0 if not int xflag,typeflag,qflag,radiusflag,rmassflag; + int replambdaHflag,moltypeHflag; int bondflag,angleflag,dihedralflag,improperflag; int nspecialflag,specialflag; int shakeflag,shakeflagflag,shakeatomflag,shaketypeflag; @@ -60,6 +61,12 @@ class Molecule : protected Pointers { double **x; // displacement of each atom from origin int *type; // type of each atom double *q; // charge on each atom + + // USER-HAdResS Package + + int *replambdaH; // replambdaH on each atom + int *moltypeH; // moltypeH on each atom + double *radius; // radius of each atom double *rmass; // mass of each atom @@ -132,6 +139,8 @@ class Molecule : protected Pointers { void coords(char *); void types(char *); void charges(char *); + void representative_atom(char *); + void moltype_atom(char *); void diameters(char *); void masses(char *); void bonds(int, char *); diff --git a/src/read_data.cpp b/src/read_data.cpp index 39afcfbcf8..e8a44fb267 100644 --- a/src/read_data.cpp +++ b/src/read_data.cpp @@ -968,7 +968,13 @@ void ReadData::header(int firstpass) } else if (strstr(line,"atom types")) { sscanf(line,"%d",&ntypes); if (addflag == NONE) atom->ntypes = ntypes + extra_atom_types; - } else if (strstr(line,"bond types")) { + } else if (strstr(line,"mol_H types")){ + + // Number of molecule types in Coarse-grain regions, USER-HADRESS Package + + sscanf(line,"%d",&atom->nmoltypesH); + } + else if (strstr(line,"bond types")) { sscanf(line,"%d",&nbondtypes); if (addflag == NONE) atom->nbondtypes = nbondtypes + extra_bond_types; } else if (strstr(line,"angle types")) {