jianmu
/
lammps
forked from lijiext/lammps
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

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

This commit is contained in:
sjplimp 2010-07-16 00:14:10 +00:00
parent eae2292f8f
commit 6f023fe957
2 changed files with 80 additions and 31 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

87
src/fix_thermal_conductivity.cpp
View File

@ -11,6 +11,11 @@
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Craig Tenney (UND) added support
for swapping atoms of different masses
------------------------------------------------------------------------- */
#include "math.h"
#include "mpi.h"
#include "string.h"
@ -225,7 +230,8 @@ void FixThermalConductivity::end_of_step()
// exchange kinetic energy between the 2 particles
// if I own both particles just swap, else point2point comm of velocities
double sbuf[3],rbuf[3];
double sbuf[4],rbuf[4],vcm[3];
double eswap = 0.0;
mine[0].proc = mine[1].proc = me;
int ilo = 0;
@ -239,45 +245,80 @@ void FixThermalConductivity::end_of_step()
MPI_Allreduce(mine,all,2,MPI_DOUBLE_INT,MPI_MINLOC,world);
if (all[0].value == BIG || all[1].value == BIG) continue;
all[0].value = -all[0].value;
e_exchange += force->mvv2e * (all[0].value - all[1].value);
if (me == all[0].proc && me == all[1].proc) {
i = index_lo[ilo++];
j = index_hi[ihi++];
sbuf[0] = v[j][0];
sbuf[1] = v[j][1];
sbuf[2] = v[j][2];
if (rmass) sbuf[3] = rmass[j];
else sbuf[3] = mass[type[j]];
rbuf[0] = v[i][0];
rbuf[1] = v[i][1];
rbuf[2] = v[i][2];
v[i][0] = v[j][0];
v[i][1] = v[j][1];
v[i][2] = v[j][2];
v[j][0] = rbuf[0];
v[j][1] = rbuf[1];
v[j][2] = rbuf[2];
if (rmass) rbuf[3] = rmass[i];
else rbuf[3] = mass[type[i]];
vcm[0] = (sbuf[3]*sbuf[0] + rbuf[3]*rbuf[0]) / (sbuf[3] + rbuf[3]);
vcm[1] = (sbuf[3]*sbuf[1] + rbuf[3]*rbuf[1]) / (sbuf[3] + rbuf[3]);
vcm[2] = (sbuf[3]*sbuf[2] + rbuf[3]*rbuf[2]) / (sbuf[3] + rbuf[3]);
v[j][0] = 2.0 * vcm[0] - sbuf[0];
v[j][1] = 2.0 * vcm[1] - sbuf[1];
v[j][2] = 2.0 * vcm[2] - sbuf[2];
eswap += sbuf[3] * (vcm[0] * (vcm[0] - sbuf[0]) +
vcm[1] * (vcm[1] - sbuf[1]) +
vcm[2] * (vcm[2] - sbuf[2]));
v[i][0] = 2.0 * vcm[0] - rbuf[0];
v[i][1] = 2.0 * vcm[1] - rbuf[1];
v[i][2] = 2.0 * vcm[2] - rbuf[2];
eswap -= rbuf[3] * (vcm[0] * (vcm[0] - rbuf[0]) +
vcm[1] * (vcm[1] - rbuf[1]) +
vcm[2] * (vcm[2] - rbuf[2]));
} else if (me == all[0].proc) {
i = index_lo[ilo++];
sbuf[0] = v[i][0];
sbuf[1] = v[i][1];
sbuf[2] = v[i][2];
MPI_Sendrecv(sbuf,3,MPI_DOUBLE,all[1].proc,0,
rbuf,3,MPI_DOUBLE,all[1].proc,0,world,&status);
v[i][0] = rbuf[0];
v[i][1] = rbuf[1];
v[i][2] = rbuf[2];
j = index_lo[ilo++];
sbuf[0] = v[j][0];
sbuf[1] = v[j][1];
sbuf[2] = v[j][2];
if (rmass) sbuf[3] = rmass[j];
else sbuf[3] = mass[type[j]];
MPI_Sendrecv(sbuf,4,MPI_DOUBLE,all[1].proc,0,
rbuf,4,MPI_DOUBLE,all[1].proc,0,world,&status);
vcm[0] = (sbuf[3]*sbuf[0] + rbuf[3]*rbuf[0]) / (sbuf[3] + rbuf[3]);
vcm[1] = (sbuf[3]*sbuf[1] + rbuf[3]*rbuf[1]) / (sbuf[3] + rbuf[3]);
vcm[2] = (sbuf[3]*sbuf[2] + rbuf[3]*rbuf[2]) / (sbuf[3] + rbuf[3]);
v[j][0] = 2.0 * vcm[0] - sbuf[0];
v[j][1] = 2.0 * vcm[1] - sbuf[1];
v[j][2] = 2.0 * vcm[2] - sbuf[2];
eswap -= sbuf[3] * (vcm[0] * (vcm[0] - sbuf[0]) +
vcm[1] * (vcm[1] - sbuf[1]) +
vcm[2] * (vcm[2] - sbuf[2]));
} else if (me == all[1].proc) {
j = index_hi[ihi++];
sbuf[0] = v[j][0];
sbuf[1] = v[j][1];
sbuf[2] = v[j][2];
MPI_Sendrecv(sbuf,3,MPI_DOUBLE,all[0].proc,0,
rbuf,3,MPI_DOUBLE,all[0].proc,0,world,&status);
v[j][0] = rbuf[0];
v[j][1] = rbuf[1];
v[j][2] = rbuf[2];
if (rmass) sbuf[3] = rmass[j];
else sbuf[3] = mass[type[j]];
MPI_Sendrecv(sbuf,4,MPI_DOUBLE,all[0].proc,0,
rbuf,4,MPI_DOUBLE,all[0].proc,0,world,&status);
vcm[0] = (sbuf[3]*sbuf[0] + rbuf[3]*rbuf[0]) / (sbuf[3] + rbuf[3]);
vcm[1] = (sbuf[3]*sbuf[1] + rbuf[3]*rbuf[1]) / (sbuf[3] + rbuf[3]);
vcm[2] = (sbuf[3]*sbuf[2] + rbuf[3]*rbuf[2]) / (sbuf[3] + rbuf[3]);
v[j][0] = 2.0 * vcm[0] - sbuf[0];
v[j][1] = 2.0 * vcm[1] - sbuf[1];
v[j][2] = 2.0 * vcm[2] - sbuf[2];
eswap += sbuf[3] * (vcm[0] * (vcm[0] - sbuf[0]) +
vcm[1] * (vcm[1] - sbuf[1]) +
vcm[2] * (vcm[2] - sbuf[2]));
}
}
// tally energy exchange from all swaps
double eswap_all;
MPI_Allreduce(&eswap,&eswap_all,1,MPI_DOUBLE,MPI_SUM,world);
e_exchange += force->mvv2e * eswap_all;
}
/* ---------------------------------------------------------------------- */

24
src/fix_viscosity.cpp
View File

@ -11,6 +11,11 @@
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Craig Tenney (UND) added support
for swapping atoms of different masses
------------------------------------------------------------------------- */
#include "math.h"
#include "mpi.h"
#include "string.h"
@ -234,7 +239,7 @@ void FixViscosity::end_of_step()
double *rmass = atom->rmass;
int ipos,ineg;
double sbuf[2],rbuf[2];
double sbuf[2],rbuf[2],vcm;
double pswap = 0.0;
mine[0].proc = mine[1].proc = me;
@ -260,9 +265,10 @@ void FixViscosity::end_of_step()
sbuf[0] = v[ineg][vdim];
if (rmass) sbuf[1] = rmass[ineg];
else sbuf[1] = mass[type[ineg]];
v[ineg][vdim] = rbuf[0] * rbuf[1]/sbuf[1];
v[ipos][vdim] = sbuf[0] * sbuf[1]/rbuf[1];
pswap += rbuf[0]*rbuf[1] - sbuf[0]*sbuf[1];
vcm = (sbuf[1]*sbuf[0] + rbuf[1]*rbuf[0]) / (sbuf[1] + rbuf[1]);
v[ineg][vdim] = 2.0 * vcm - sbuf[0];
v[ipos][vdim] = 2.0 * vcm - rbuf[0];
pswap += rbuf[1] * (vcm - rbuf[0]) - sbuf[1] * (vcm - sbuf[0]);
} else if (me == all[0].proc) {
ipos = pos_index[ipositive++];
@ -271,8 +277,9 @@ void FixViscosity::end_of_step()
else sbuf[1] = mass[type[ipos]];
MPI_Sendrecv(sbuf,2,MPI_DOUBLE,all[1].proc,0,
rbuf,2,MPI_DOUBLE,all[1].proc,0,world,&status);
v[ipos][vdim] = rbuf[0] * rbuf[1]/sbuf[1];
pswap += sbuf[0]*sbuf[1];
vcm = (sbuf[1]*sbuf[0] + rbuf[1]*rbuf[0]) / (sbuf[1] + rbuf[1]);
v[ipos][vdim] = 2.0 * vcm - sbuf[0];
pswap += sbuf[1] * (vcm - sbuf[0]);
} else if (me == all[1].proc) {
ineg = neg_index[inegative++];
@ -281,8 +288,9 @@ void FixViscosity::end_of_step()
else sbuf[1] = mass[type[ineg]];
MPI_Sendrecv(sbuf,2,MPI_DOUBLE,all[0].proc,0,
rbuf,2,MPI_DOUBLE,all[0].proc,0,world,&status);
v[ineg][vdim] = rbuf[0] * rbuf[1]/sbuf[1];
pswap -= sbuf[0]*sbuf[1];
vcm = (sbuf[1]*sbuf[0] + rbuf[1]*rbuf[0]) / (sbuf[1] + rbuf[1]);
v[ineg][vdim] = 2.0 * vcm - sbuf[0];
pswap -= sbuf[1] * (vcm - sbuf[0]);
}
}