git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@1388 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp
2008-01-18 17:18:46 +00:00
parent 89beb4487e
commit 17416f49f4
8 changed files with 36 additions and 42 deletions

View File

@ -824,7 +824,7 @@ int PairMEAM::pack_reverse_comm(int n, int first, double *buf)
buf[m++] = t_ave[i][2]; buf[m++] = t_ave[i][2];
} }
return m; return 27;
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

View File

@ -365,8 +365,8 @@ void PairBuckCoul::write_restart(FILE *fp)
fwrite(&setflag[i][j],sizeof(int),1,fp); fwrite(&setflag[i][j],sizeof(int),1,fp);
if (setflag[i][j]) { if (setflag[i][j]) {
fwrite(&buck_a_read[i][j],sizeof(double),1,fp); fwrite(&buck_a_read[i][j],sizeof(double),1,fp);
fwrite(&buck_c_read[i][j],sizeof(double),1,fp);
fwrite(&buck_rho_read[i][j],sizeof(double),1,fp); fwrite(&buck_rho_read[i][j],sizeof(double),1,fp);
fwrite(&buck_c_read[i][j],sizeof(double),1,fp);
fwrite(&cut_buck_read[i][j],sizeof(double),1,fp); fwrite(&cut_buck_read[i][j],sizeof(double),1,fp);
} }
} }
@ -391,13 +391,13 @@ void PairBuckCoul::read_restart(FILE *fp)
if (setflag[i][j]) { if (setflag[i][j]) {
if (me == 0) { if (me == 0) {
fread(&buck_a_read[i][j],sizeof(double),1,fp); fread(&buck_a_read[i][j],sizeof(double),1,fp);
fread(&buck_c_read[i][j],sizeof(double),1,fp);
fread(&buck_rho_read[i][j],sizeof(double),1,fp); fread(&buck_rho_read[i][j],sizeof(double),1,fp);
fread(&buck_c_read[i][j],sizeof(double),1,fp);
fread(&cut_buck_read[i][j],sizeof(double),1,fp); fread(&cut_buck_read[i][j],sizeof(double),1,fp);
} }
MPI_Bcast(&buck_a_read[i][j],1,MPI_DOUBLE,0,world); MPI_Bcast(&buck_a_read[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&buck_c_read[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&buck_rho_read[i][j],1,MPI_DOUBLE,0,world); MPI_Bcast(&buck_rho_read[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&buck_c_read[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&cut_buck_read[i][j],1,MPI_DOUBLE,0,world); MPI_Bcast(&cut_buck_read[i][j],1,MPI_DOUBLE,0,world);
} }
} }

View File

@ -167,6 +167,8 @@ void DisplaceAtoms::command(int narg, char **arg)
} }
// move atoms randomly // move atoms randomly
// makes atom result independent of what proc owns it via random->reset()
if (style == RANDOM) { if (style == RANDOM) {
RanPark *random = new RanPark(lmp,1); RanPark *random = new RanPark(lmp,1);

View File

@ -95,54 +95,44 @@ void RanPark::reset(int seed_init)
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
reset the seed based on atom position within box and ibase = caller seed reset the seed based on atom coords and ibase = caller seed
combine 3 RNGs based on fractional position in box into one new seed use hash function, treating user seed and coords as sequence of input ints
this is Jenkins One-at-a-time hash, see Wikipedia entry on hash tables
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void RanPark::reset(int ibase, double *coord) void RanPark::reset(int ibase, double *coord)
{ {
// for orthogonal box, lamda = fractional position in box int i;
// for triclinic box, convert to lamda coords
double lamda[3]; char *str = (char *) &ibase;
int n = sizeof(int);
if (domain->triclinic == 0) { unsigned int hash = 0;
lamda[0] = (coord[0] - domain->boxlo[0]) / domain->prd[0]; for (i = 0; i < n; i++) {
lamda[1] = (coord[1] - domain->boxlo[1]) / domain->prd[1]; hash += str[i];
lamda[2] = (coord[2] - domain->boxlo[2]) / domain->prd[2]; hash += (hash << 10);
} else domain->x2lamda(coord,lamda); hash ^= (hash >> 6);
}
// seed 1,2,3 = combination of atom coord in each dim and user-input seed str = (char *) coord;
// map geometric extent into range of each of 3 RNGs n = 3 * sizeof(double);
// warm-up each RNG by calling it twice for (i = 0; i < n; i++) {
hash += str[i];
hash += (hash << 10);
hash ^= (hash >> 6);
}
int seed1,seed2,seed3; hash += (hash << 3);
hash ^= (hash >> 11);
hash += (hash << 15);
seed1 = static_cast<int> (lamda[0] * IM1); // keep 31 bits of unsigned int as new seed
seed1 = (seed1+ibase) % IM1;
seed1 = (seed1*IA1+IC1) % IM1;
seed1 = (seed1*IA1+IC1) % IM1;
seed2 = static_cast<int> (lamda[1] * IM2); seed = hash & 0x7ffffff;
seed2 = (seed2+ibase) % IM2;
seed2 = (seed2*IA2+IC2) % IM2;
seed2 = (seed2*IA2+IC2) % IM2;
seed3 = static_cast<int> (lamda[2] * IM3); // warm up the RNG
seed3 = (seed3+ibase) % IM3;
seed3 = (seed3*IA3+IC3) % IM3;
seed3 = (seed3*IA3+IC3) % IM3;
// fraction = 0-1 with giving each dim an equal weighting for (i = 0; i < 5; i++) uniform();
// use fraction to reset Park/Miller RNG seed
// warm-up master RNG with new seed by calling it twice
double fraction = 1.0*seed1/(3*IM1) + 1.0*seed2/(3*IM2) + 1.0*seed3/(3*IM3);
seed = static_cast<int> (fraction*IM) + 1;
if (seed >= IM) seed = IM-1;
uniform();
uniform();
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

View File

@ -19,6 +19,7 @@
namespace LAMMPS_NS { namespace LAMMPS_NS {
class RanPark : protected Pointers { class RanPark : protected Pointers {
friend class Set;
public: public:
RanPark(class LAMMPS *, int); RanPark(class LAMMPS *, int);
double uniform(); double uniform();

View File

@ -324,7 +324,7 @@ void Set::set(int keyword)
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
set an owned atom property randomly set an owned atom property randomly
set seed based on atom tag set seed based on atom tag
makes atom's result independent of what proc owns it make atom result independent of what proc owns it
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Set::setrandom(int keyword) void Set::setrandom(int keyword)

View File

@ -170,6 +170,7 @@ void Velocity::create(int narg, char **arg)
// will never produce same V, independent of P // will never produce same V, independent of P
// GEOM = only loop over my atoms // GEOM = only loop over my atoms
// choose RNG for each atom based on its xyz coord (geometry) // choose RNG for each atom based on its xyz coord (geometry)
// via random->reset()
// will always produce same V, independent of P // will always produce same V, independent of P
// adjust by factor for atom mass // adjust by factor for atom mass
// for 2d, set Vz to 0.0 // for 2d, set Vz to 0.0

View File

@ -100,7 +100,7 @@ void Verlet::setup()
ev_set(update->ntimestep); ev_set(update->ntimestep);
force_clear(); force_clear();
if (force->pair) force->pair->compute(eflag,vflag); if (force->pair) force->pair->compute(eflag,vflag);
if (atom->molecular) { if (atom->molecular) {