restore bugfixes and updates that were lost. flag time dependet. correct use of citeme. (b27179bb) · Commits · 郑智淋 / lammps

examples/USER/misc/ti/in.ti_spring

100755 → 100644

+0 −0

File mode changed from 100755 to 100644.

src/USER-MISC/fix_ti_spring.cpp

+51 −43

Original line number	Diff line number	Diff line
		@@ -18,8 +18,8 @@
		Maurice de Koning (Unicamp/Brazil) - dekoning@ifi.unicamp.br
		------------------------------------------------------------------------- */

		#include "stdlib.h"
		#include "string.h"
		#include <stdlib.h>
		#include <string.h>
		#include "fix_ti_spring.h"
		#include "atom.h"
		#include "update.h"
		@@ -27,12 +27,13 @@
		#include "respa.h"
		#include "memory.h"
		#include "error.h"
		#include "citeme.h"
		#include "force.h"

		using namespace LAMMPS_NS;
		using namespace FixConst;

		static const char cite_ti_spring[] =
		static const char cite_fix_ti_spring[] =
		"ti/spring command:\n\n"
		"@article{freitas2016,\n"
		" author={Freitas, Rodrigo and Asta, Mark and de Koning, Maurice},\n"
		@@ -49,6 +50,8 @@ static const char cite_ti_spring[] =
		FixTISpring::FixTISpring(LAMMPS lmp, int narg, char *arg) :
		Fix(lmp, narg, arg)
		{
		if (lmp->citeme) lmp->citeme->add(cite_fix_ti_spring);

		if (narg < 6 \|\| narg > 8)
		error->all(FLERR,"Illegal fix ti/spring command");

		@@ -62,22 +65,25 @@ FixTISpring::FixTISpring(LAMMPS lmp, int narg, char *arg) :
		extscalar = 1;
		extvector = 1;

		// disallow resetting the time step, while this fix is defined
		time_depend = 1;

		// Spring constant.
		k = force->numeric(FLERR,arg[3]);
		if (k <= 0.0) error->all(FLERR,"Illegal fix ti/spring command");

		// Perform initial allocation of atom-based array.
		// Registar with Atom class.
		// Perform initial allocation of atom-based array
		// Register with Atom class
		xoriginal = NULL;
		grow_arrays(atom->nmax);
		atom->add_callback(0);
		atom->add_callback(1);

		// xoriginal = initial unwrapped positions of atom.
		// xoriginal = initial unwrapped positions of atoms

		double **x = atom->x;
		int *mask = atom->mask;
		tagint *image = atom->image;
		imageint *image = atom->image;
		int nlocal = atom->nlocal;

		for (int i = 0; i < nlocal; i++) {
		@@ -86,16 +92,16 @@ FixTISpring::FixTISpring(LAMMPS lmp, int narg, char *arg) :
		}

		// Time variables.
		t_switch = atoi(arg[4]); // Switching time.
		t_equil = atoi(arg[5]); // Equilibration time.
		t0 = update->ntimestep; // Initial time.
		if (t_switch < 0.0) error->all(FLERR,"Illegal fix ti/spring command");
		if (t_equil < 0.0) error->all(FLERR,"Illegal fix ti/spring command");
		t0 = update->ntimestep; // timestep of original/starting coordinates
		t_switch = force->bnumeric(FLERR,arg[4]); // Number of steps for switching
		t_equil = force->bnumeric(FLERR,arg[5]); // Number of steps for equilibration
		if ((t_switch <= 0) \|\| (t_equil < 0))
		error->all(FLERR,"Illegal fix ti/spring command");

		// Coupling parameter initialization.
		// Coupling parameter initialization
		sf = 1;
		if (narg > 6) {
		if (strcmp(arg[6], "function") == 0) sf = atoi(arg[7]);
		if (strcmp(arg[6], "function") == 0) sf = force->inumeric(FLERR,arg[7]);
		else error->all(FLERR,"Illegal fix ti/spring switching function");
		if ((sf!=1) && (sf!=2))
		error->all(FLERR,"Illegal fix ti/spring switching function");
		@@ -163,14 +169,13 @@ void FixTISpring::min_setup(int vflag)

		void FixTISpring::post_force(int vflag)
		{
		// If on the first equilibration do not calculate forces.
		int t = update->ntimestep - t0;
		if(t < t_equil) return;
		// do not calculate forces during equilibration
		if ((update->ntimestep - t0) < t_equil) return;

		double **x = atom->x;
		double **f = atom->f;
		int *mask = atom->mask;
		tagint *image = atom->image;
		imageint *image = atom->image;
		int nlocal = atom->nlocal;

		double dx, dy, dz;
		@@ -211,17 +216,20 @@ void FixTISpring::min_post_force(int vflag)

		void FixTISpring::initial_integrate(int vflag)
		{
		// Update the coupling parameter value.
		double t = update->ntimestep - (t0+t_equil);
		// Update the coupling parameter value if needed
		if ((update->ntimestep - t0) < t_equil) return;

		const bigint t = update->ntimestep - (t0+t_equil);
		const double r_switch = 1.0/t_switch;

		if ( (t >= 0) && (t <= t_switch) ) {
		lambda = switch_func(t/t_switch);
		dlambda = dswitch_func(t/t_switch);
		lambda = switch_func(t*r_switch);
		dlambda = dswitch_func(t*r_switch);
		}

		if ( (t >= t_equil+t_switch) && (t <= (t_equil+2*t_switch)) ) {
		lambda = switch_func(1.0 - (t - t_switch - t_equil)/t_switch );
		dlambda = - dswitch_func(1.0 - (t - t_switch - t_equil)/t_switch );
		lambda = switch_func(1.0 - (t - t_switch - t_equil)*r_switch);
		dlambda = - dswitch_func(1.0 - (t - t_switch - t_equil)*r_switch);
		}
		}

		@@ -352,7 +360,7 @@ int FixTISpring::size_restart(int nlocal)
		}

		/* ----------------------------------------------------------------------
		Switching function.
		Switching function
		------------------------------------------------------------------------- */

		double FixTISpring::switch_func(double t)
		@@ -365,7 +373,7 @@ double FixTISpring::switch_func(double t)
		}

		/* ----------------------------------------------------------------------
		Switching function derivative.
		Switching function derivative
		------------------------------------------------------------------------- */

		double FixTISpring::dswitch_func(double t)

src/USER-MISC/fix_ti_spring.h

+3 −3

Original line number	Diff line number	Diff line
		@@ -66,9 +66,9 @@ class FixTISpring : public Fix {
		double lambda; // Coupling parameter.
		double dlambda; // Lambda variation with t.
		double linfo[2]; // Current lambda status.
		int t_switch; // Total switching steps.
		int t_equil; // Equilibration time.
		int t0; // Initial time.
		bigint t_switch; // Total switching steps.
		bigint t_equil; // Equilibration time.
		bigint t0; // Initial time.
		int sf; // Switching function option.
		int nlevels_respa;
		};

Admin message