Merge pull request #1628 from erozic/feature-wall-region-morse (82b3fad1) · Commits · 郑智淋 / lammps

doc/src/fix_wall_region.txt

+31 −14

Original line number	Diff line number	Diff line
		@@ -10,19 +10,27 @@ fix wall/region command :h3

		[Syntax:]

		fix ID group-ID wall/region region-ID style epsilon sigma cutoff :pre
		fix ID group-ID wall/region region-ID style args ... cutoff :pre

		ID, group-ID are documented in "fix"_fix.html command
		wall/region = style name of this fix command
		region-ID = region whose boundary will act as wall
		style = {lj93} or {lj126} or {lj1043} or {colloid} or {harmonic}
		ID, group-ID are documented in "fix"_fix.html command :ulb,l
		wall/region = style name of this fix command :l
		region-ID = region whose boundary will act as wall :l
		style = {lj93} or {lj126} or {lj1043} or {colloid} or {harmonic} or {morse} :l
		args for styles {lj93} or {lj126} or {lj1043} or {colloid} or {harmonic} = :l
		epsilon = strength factor for wall-particle interaction (energy or energy/distance^2 units)
		sigma = size factor for wall-particle interaction (distance units)
		cutoff = distance from wall at which wall-particle interaction is cut off (distance units) :ul
		sigma = size factor for wall-particle interaction (distance units) :pre
		args for style {morse} = :l
		D_0 = depth of the potential (energy units)
		alpha = width parameter (1/distance units)
		r_0 = distance of the potential minimum from wall position (distance units) :pre
		cutoff = distance from wall at which wall-particle interaction is cut off (distance units) :l
		:ule

		[Examples:]

		fix wall all wall/region mySphere lj93 1.0 1.0 2.5 :pre
		fix wall all wall/region mySphere lj93 1.0 1.0 2.5
		fix wall all wall/region mySphere harmonic 1.0 0.0 2.5
		fix wall all wall/region box_top morse 1.0 1.0 1.5 3.0 :pre

		[Description:]

		@@ -122,15 +130,22 @@ the "pair_style colloid"_pair_colloid.html potential:
		:c,image(Eqs/fix_wall_colloid.jpg)

		For style {wall/harmonic}, the energy E is given by a harmonic spring
		potential:
		potential (the distance parameter is ignored):

		:c,image(Eqs/fix_wall_harmonic.jpg)

		For style {wall/morse}, the energy E is given by the Morse potential:

		:c,image(Eqs/pair_morse.jpg)

		Unlike other styles, this requires three parameters ({D_0}, {alpha}, {r_0}
		in this order) instead of two like for the other wall styles.

		In all cases, {r} is the distance from the particle to the region
		surface, and Rc is the {cutoff} distance at which the particle and
		surface no longer interact. The energy of the wall potential is
		shifted so that the wall-particle interaction energy is 0.0 at the
		cutoff distance.
		surface no longer interact. The cutoff is always the last argument.
		The energy of the wall potential is shifted so that the wall-particle
		interaction energy is 0.0 at the cutoff distance.

		For a full description of these wall styles, see fix_style
		"wall"_fix_wall.html
		@@ -179,7 +194,9 @@ option for this fix.

		"fix wall/lj93"_fix_wall.html,
		"fix wall/lj126"_fix_wall.html,
		"fix wall/lj1043"_fix_wall.html,
		"fix wall/colloid"_fix_wall.html,
		"fix wall/harmonic"_fix_wall.html,
		"fix wall/gran"_fix_wall_gran.html

		[Default:] none

src/fix_wall_region.cpp

+39 −7

Original line number	Diff line number	Diff line
		@@ -28,7 +28,7 @@ using namespace LAMMPS_NS;
		using namespace FixConst;
		using namespace MathConst;

		enum{LJ93,LJ126,LJ1043,COLLOID,HARMONIC};
		enum{LJ93,LJ126,LJ1043,COLLOID,HARMONIC,MORSE};

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

		@@ -36,7 +36,7 @@ FixWallRegion::FixWallRegion(LAMMPS lmp, int narg, char *arg) :
		Fix(lmp, narg, arg),
		idregion(NULL)
		{
		if (narg != 8) error->all(FLERR,"Illegal fix wall/region command");
		if (narg < 8) error->all(FLERR,"Illegal fix wall/region command");

		scalar_flag = 1;
		vector_flag = 1;
		@@ -62,13 +62,28 @@ FixWallRegion::FixWallRegion(LAMMPS lmp, int narg, char *arg) :
		else if (strcmp(arg[4],"lj1043") == 0) style = LJ1043;
		else if (strcmp(arg[4],"colloid") == 0) style = COLLOID;
		else if (strcmp(arg[4],"harmonic") == 0) style = HARMONIC;
		else if (strcmp(arg[4],"morse") == 0) style = MORSE;
		else error->all(FLERR,"Illegal fix wall/region command");

		if (style != COLLOID) dynamic_group_allow = 1;

		if (style == MORSE) {
		if (narg != 9)
		error->all(FLERR,"Illegal fix wall/region command");

		epsilon = force->numeric(FLERR,arg[5]);
		alpha = force->numeric(FLERR,arg[6]);
		sigma = force->numeric(FLERR,arg[7]);
		cutoff = force->numeric(FLERR,arg[8]);

		} else {
		if (narg != 8)
		error->all(FLERR,"Illegal fix wall/region command");

		epsilon = force->numeric(FLERR,arg[5]);
		sigma = force->numeric(FLERR,arg[6]);
		cutoff = force->numeric(FLERR,arg[7]);
		}

		if (cutoff <= 0.0) error->all(FLERR,"Fix wall/region cutoff <= 0.0");

		@@ -154,12 +169,15 @@ void FixWallRegion::init()
		coeff5 = coeff1 * 10.0;
		coeff6 = coeff2 * 4.0;
		coeff7 = coeff3 * 3.0;

		double rinv = 1.0/cutoff;
		double r2inv = rinv*rinv;
		double r4inv = r2inv*r2inv;
		offset = coeff1r4invr4invr2inv - coeff2r4inv -
		coeff3*pow(cutoff+coeff4,-3.0);
		} else if (style == MORSE) {
		coeff1 = 2 * epsilon * alpha;
		double alpha_dr = -alpha * (cutoff - sigma);
		offset = epsilon * (exp(2.0alpha_dr) - 2.0exp(alpha_dr));
		} else if (style == COLLOID) {
		coeff1 = -4.0/315.0 * epsilon * pow(sigma,6.0);
		coeff2 = -2.0/3.0 * epsilon;
		@@ -250,6 +268,7 @@ void FixWallRegion::post_force(int vflag)
		if (style == LJ93) lj93(region->contact[m].r);
		else if (style == LJ126) lj126(region->contact[m].r);
		else if (style == LJ1043) lj1043(region->contact[m].r);
		else if (style == MORSE) morse(region->contact[m].r);
		else if (style == COLLOID) colloid(region->contact[m].r,radius[i]);
		else harmonic(region->contact[m].r);

		@@ -372,6 +391,19 @@ void FixWallRegion::lj1043(double r)
		coeff3*pow(r+coeff4,-3.0) - offset;
		}

		/* ----------------------------------------------------------------------
		Morse interaction for particle with wall
		compute eng and fwall = magnitude of wall force
		------------------------------------------------------------------------- */

		void FixWallRegion::morse(double r)
		{
		double dr = r - sigma;
		double dexp = exp(-alpha * dr);
		fwall = coeff1 * (dexp*dexp - dexp) / r;
		eng = epsilon * (dexpdexp - 2.0dexp) - offset;
		}

		/* ----------------------------------------------------------------------
		colloid interaction for finite-size particle of rad with wall
		compute eng and fwall = magnitude of wall force

src/fix_wall_region.h

+2 −0

Original line number	Diff line number	Diff line
		@@ -41,6 +41,7 @@ class FixWallRegion : public Fix {
		private:
		int style,iregion;
		double epsilon,sigma,cutoff;
		double alpha;
		int eflag;
		double ewall[4],ewall_all[4];
		int ilevel_respa;
		@@ -53,6 +54,7 @@ class FixWallRegion : public Fix {
		void lj93(double);
		void lj126(double);
		void lj1043(double);
		void morse(double);
		void colloid(double, double);
		void harmonic(double);
		};

Admin message