add virial support to fix addforce (8776b810) · Commits · 郑智淋 / lammps

src/fix_addforce.cpp

+38 −6

Original line number	Diff line number	Diff line
		@@ -51,6 +51,7 @@ FixAddForce::FixAddForce(LAMMPS lmp, int narg, char *arg) :
		extvector = 1;
		respa_level_support = 1;
		ilevel_respa = 0;
		virial_flag = 1;

		xstr = ystr = zstr = NULL;

		@@ -237,10 +238,16 @@ void FixAddForce::post_force(int vflag)
		double **f = atom->f;
		int *mask = atom->mask;
		imageint *image = atom->image;
		double v[6];
		int nlocal = atom->nlocal;

		if (update->ntimestep % nevery) return;

		// energy and virial setup

		if (vflag) v_setup(vflag);
		else evflag = 0;

		if (lmp->kokkos)
		atom->sync_modify(Host, (unsigned int) (F_MASK \| MASK_MASK),
		(unsigned int) F_MASK);
		@@ -283,6 +290,14 @@ void FixAddForce::post_force(int vflag)
		f[i][0] += xvalue;
		f[i][1] += yvalue;
		f[i][2] += zvalue;
		if (evflag) {
		v[0] = xvalue * unwrap[0];
		v[1] = yvalue * unwrap[1];
		v[2] = zvalue * unwrap[2];
		v[3] = xvalue * unwrap[1];
		v[4] = xvalue * unwrap[2];
		v[5] = yvalue * unwrap[2];
		v_tally(i,v);
		}

		// variable force, wrap with clear/add
		@@ -290,6 +305,7 @@ void FixAddForce::post_force(int vflag)
		// wrap with clear/add

		} else {
		double unwrap[3];

		modify->clearstep_compute();

		@@ -310,16 +326,32 @@ void FixAddForce::post_force(int vflag)
		for (int i = 0; i < nlocal; i++)
		if (mask[i] & groupbit) {
		if (region && !region->match(x[i][0],x[i][1],x[i][2])) continue;
		domain->unmap(x[i],image[i],unwrap);
		if (xstyle == ATOM) xvalue = sforce[i][0];
		if (ystyle == ATOM) yvalue = sforce[i][1];
		if (zstyle == ATOM) zvalue = sforce[i][2];

		if (estyle == ATOM) foriginal[0] += sforce[i][3];
		else {
		if (xstyle) foriginal[0] -= xvalue*unwrap[0];
		if (ystyle) foriginal[0] -= yvalue*unwrap[1];
		if (zstyle) foriginal[0] -= zvalue*unwrap[2];
		}
		foriginal[1] += f[i][0];
		foriginal[2] += f[i][1];
		foriginal[3] += f[i][2];
		if (xstyle == ATOM) f[i][0] += sforce[i][0];
		else if (xstyle) f[i][0] += xvalue;
		if (ystyle == ATOM) f[i][1] += sforce[i][1];
		else if (ystyle) f[i][1] += yvalue;
		if (zstyle == ATOM) f[i][2] += sforce[i][2];
		else if (zstyle) f[i][2] += zvalue;
		if (xstyle) f[i][0] += xvalue;
		if (ystyle) f[i][1] += yvalue;
		if (zstyle) f[i][2] += zvalue;
		if (evflag) {
		v[0] = xstyle ? fx*unwrap[0] : 0.0;
		v[1] = ystyle ? fy*unwrap[1] : 0.0;
		v[2] = zstyle ? fz*unwrap[2] : 0.0;
		v[3] = xstyle ? fx*unwrap[1] : 0.0;
		v[4] = xstyle ? fx*unwrap[2] : 0.0;
		v[5] = ystyle ? fy*unwrap[2] : 0.0;
		v_tally(i, v);
		}
		}
		}
		}

Admin message