Stopped working on gaussian bump. (7adc7f02) · Commits · 郑智淋 / lammps

src/USER-MANIFOLD/fix_nve_manifold_rattle.cpp

+7 −0

Original line number	Diff line number	Diff line
		@@ -115,6 +115,13 @@ FixNVEManifoldRattle::FixNVEManifoldRattle( LAMMPS lmp, int &narg, char *arg,
		error->all(FLERR, "Error creating manifold arg arrays");
		}

		// Check if you have enough args:
		if( 6 + nvars > narg ){
		char msg[2048];
		sprintf(msg, "Not enough args for manifold %s, %d expected but got %d\n",
		ptr_m->id(), nvars, narg - 6);
		error->all(FLERR, msg);
		}
		// Loop over manifold args:
		for( int i = 0; i < nvars; ++i ){
		int len = 0, offset = 0;

src/USER-MANIFOLD/manifold_factory.cpp

+3 −1

Original line number	Diff line number	Diff line
		@@ -37,6 +37,7 @@
		#include "manifold_cylinder_dent.h"
		#include "manifold_dumbbell.h"
		#include "manifold_ellipsoid.h"
		#include "manifold_gaussian_bump.h"
		#include "manifold_plane.h"
		#include "manifold_plane_wiggle.h"
		#include "manifold_sphere.h"
		@@ -58,6 +59,7 @@ manifold* LAMMPS_NS::user_manifold::create_manifold(const char *mname,
		make_manifold_if<manifold_cylinder_dent> ( &man, mname, lmp, narg, arg );
		make_manifold_if<manifold_dumbbell> ( &man, mname, lmp, narg, arg );
		make_manifold_if<manifold_ellipsoid> ( &man, mname, lmp, narg, arg );
		make_manifold_if<manifold_gaussian_bump> ( &man, mname, lmp, narg, arg );
		make_manifold_if<manifold_plane> ( &man, mname, lmp, narg, arg );
		make_manifold_if<manifold_plane_wiggle> ( &man, mname, lmp, narg, arg );
		make_manifold_if<manifold_sphere> ( &man, mname, lmp, narg, arg );

src/USER-MANIFOLD/manifold_gaussian_bump.cpp

0 → 100644

+138 −0

Original line number	Diff line number	Diff line
		#include "manifold_gaussian_bump.h"

		using namespace LAMMPS_NS;
		using namespace user_manifold;

		// The constraint is z = f(r = sqrt( x^2 + y^2) )
		// f(x) = A*exp( -x^2 / 2 l^2 ) if x < rc1
		// = a + bx + cx*2 + dx**3 if rc1 <= x < rc2
		// = 0 if x >= rc2
		//
		double manifold_gaussian_bump::g( const double *x )
		{
		double xf[3];
		xf[0] = x[0];
		xf[1] = x[1];
		xf[2] = 0.0;

		double x2 = dot(xf,xf);
		if( x2 < rc12 ){
		return x[2] - gaussian_bump_x2( x2 );
		}else if( x2 < rc22 ){
		double rr = sqrt(x2);
		double xi = rr - rc1;
		xi *= inv_dr;
		double xi2 = x2 * inv_dr*inv_dr;
		double xi3 = xi*xi2;
		return x[2] - ( aa + bbxi + ccxi2 + dd*xi3 );

		}else{
		return x[2];
		}
		}

		void manifold_gaussian_bump::n( const double x, double nn )
		{
		double xf[3];
		xf[0] = x[0];
		xf[1] = x[1];
		xf[2] = 0.0;
		nn[2] = 1.0;

		double x2 = dot(xf,xf);

		if( x2 < rc12 ){
		double factor = gaussian_bump_x2(x2);
		factor /= (ll*ll);
		nn[0] = factor * x[0];
		nn[1] = factor * x[1];
		}else if( x2 < rc22 ){
		double rr = sqrt(x2);
		double xi = rr - rc1;
		xi *= inv_dr;
		double xi2 = x2 * inv_dr*inv_dr;
		double der = bb + 2ccxi + 3ddxi2;

		nn[0] = -der * x[0] / rr;
		nn[1] = -der * x[1] / rr;
		}else{
		nn[0] = nn[1] = 0.0;
		}
		}

		double manifold_gaussian_bump::g_and_n( const double x, double nn )
		{
		double xf[3];
		xf[0] = x[0];
		xf[1] = x[1];
		xf[2] = 0.0;
		nn[2] = 1.0;

		double x2 = dot(xf,xf);
		if( x2 < rc12 ){
		double gb = gaussian_bump_x2(x2);
		double factor = gb / (ll*ll);
		nn[0] = factor * x[0];
		nn[1] = factor * x[1];

		return x[2] - gb;
		}else if( x2 < rc22 ){

		double rr = sqrt(x2);
		double xi = rr - rc1;
		xi *= inv_dr;
		double xi2 = x2 * inv_dr*inv_dr;
		double xi3 = xi*xi2;

		double der = bb + 2ccxi + 3ddxi2;

		nn[0] = -der * x[0] / rr;
		nn[1] = -der * x[1] / rr;


		return x[2] - ( aa + bbxi + ccxi2 + dd*xi3 );

		}else{
		nn[0] = nn[1] = 0.0;
		return x[2];
		}

		}


		void manifold_gaussian_bump::post_param_init()
		{
		// Read in the params:
		AA = params[0];
		ll = params[1];
		rc1 = params[2];
		rc2 = params[3];

		ll2 = 2.0llll;

		f_at_rc = gaussian_bump_x2 ( rc12 );
		fp_at_rc = gaussian_bump_der( rc12 );

		rc12 = rc1*rc1;
		rc22 = rc2*rc2;
		dr = rc2 - rc1;
		inv_dr = 1.0 / dr;
		}


		double manifold_gaussian_bump::gaussian_bump( double x )
		{
		double x2 = x*x;
		return gaussian_bump_x2( x2 );
		}

		double manifold_gaussian_bump::gaussian_bump_x2( double x2 )
		{
		return AA*exp( -x2 / ll2 );
		}

		double manifold_gaussian_bump::gaussian_bump_der( double x )
		{
		double x2 = x*x;
		return gaussian_bump_x2( x2 )( -x/(llll) );
		}

src/USER-MANIFOLD/manifold_gaussian_bump.h

0 → 100644

+79 −0

Original line number	Diff line number	Diff line
		/* ----------------------------------------------------------------------
		LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
		http://lammps.sandia.gov, Sandia National Laboratories
		Steve Plimpton, sjplimp@sandia.gov

		Copyright (2003) Sandia Corporation. Under the terms of Contract
		DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
		certain rights in this software. This software is distributed under
		the GNU General Public License.

		See the README file in the top-level LAMMPS directory.
		-----------------------------------------------------------------------

		This file is a part of the USER-MANIFOLD package.

		This package allows LAMMPS to perform MD simulations of particles
		constrained on a manifold (i.e., a 2D subspace of the 3D simulation
		box). It achieves this using the RATTLE constraint algorithm applied
		to single-particle constraint functions g(xi,yi,zi) = 0 and their
		derivative (i.e. the normal of the manifold) n = grad(g).

		It is very easy to add your own manifolds to the current zoo
		(we now have sphere, a dendritic spine approximation, a 2D plane (for
		testing purposes) and a wave-y plane.
		See the README file for more info.

		Stefan Paquay, s.paquay@tue.nl
		Applied Physics/Theory of Polymers and Soft Matter,
		Eindhoven University of Technology (TU/e), The Netherlands

		Thanks to Remy Kusters at TU/e for testing.

		This software is distributed under the GNU General Public License.

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

		#ifndef LMP_MANIFOLD_GAUSSIAN_BUMP_H
		#define LMP_MANIFOLD_GAUSSIAN_BUMP_H

		#include "manifold.h"

		namespace LAMMPS_NS {

		namespace user_manifold {

		// A Gaussian bump with a smoothed decay to flat 2D.
		class manifold_gaussian_bump : public manifold {
		public:
		enum { NPARAMS = 4 };
		manifold_gaussian_bump(class LAMMPS, int, char *) : manifold(lmp) {}
		virtual ~manifold_gaussian_bump(){}

		virtual double g( const double * );
		virtual void n( const double , double );

		// Variant of g that computes n at the same time.
		virtual double g_and_n( const double x, double nn );

		static const char* type(){ return "gaussian_bump"; }
		virtual const char *id() { return type(); }

		virtual int nparams(){ return NPARAMS; }
		virtual void post_param_init();
		private:
		// Some private constants:
		double aa, bb, cc, dd, AA, ll, ll2, f_at_rc, fp_at_rc;
		double rc1, rc2, rc12, rc22, dr, inv_dr;

		double gaussian_bump ( double );
		double gaussian_bump_x2 ( double );
		double gaussian_bump_der( double );

		};
		}

		}


		#endif // LMP_MANIFOLD_GAUSSIAN_BUMP_H

Admin message