(35e73534) · Commits · 郑智淋 / lammps

src/pair_coul_cut.h

+7 −7

Original line number	Diff line number	Diff line
		@@ -21,19 +21,19 @@ namespace LAMMPS_NS {
		class PairCoulCut : public Pair {
		public:
		PairCoulCut(class LAMMPS *);
		~PairCoulCut();
		void compute(int, int);
		void settings(int, char **);
		virtual ~PairCoulCut();
		virtual void compute(int, int);
		virtual void settings(int, char **);
		void coeff(int, char **);
		void init_style();
		double init_one(int, int);
		void write_restart(FILE *);
		void read_restart(FILE *);
		void write_restart_settings(FILE *);
		void read_restart_settings(FILE *);
		double single(int, int, int, int, double, double, double, double &);
		virtual void write_restart_settings(FILE *);
		virtual void read_restart_settings(FILE *);
		virtual double single(int, int, int, int, double, double, double, double &);

		private:
		protected:
		double cut_global;
		double **cut;

src/pair_coul_debye.cpp

0 → 100644

+200 −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.
		------------------------------------------------------------------------- */

		#include "math.h"
		#include "stdio.h"
		#include "stdlib.h"
		#include "string.h"
		#include "pair_coul_debye.h"
		#include "atom.h"
		#include "comm.h"
		#include "force.h"
		#include "neighbor.h"
		#include "neigh_list.h"
		#include "memory.h"
		#include "error.h"

		using namespace LAMMPS_NS;

		#define MIN(a,b) ((a) < (b) ? (a) : (b))
		#define MAX(a,b) ((a) > (b) ? (a) : (b))

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

		PairCoulDebye::PairCoulDebye(LAMMPS *lmp) : PairCoulCut(lmp) {}

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

		PairCoulDebye::~PairCoulDebye()
		{
		if (allocated) {
		memory->destroy_2d_int_array(setflag);
		memory->destroy_2d_double_array(cutsq);

		memory->destroy_2d_double_array(cut);
		}
		}

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

		void PairCoulDebye::compute(int eflag, int vflag)
		{
		int i,j,ii,jj,inum,jnum,itype,jtype;
		double qtmp,xtmp,ytmp,ztmp,delx,dely,delz,ecoul,fpair;
		double rsq,r2inv,r,rinv,forcecoul,factor_coul,screening;
		int ilist,jlist,numneigh,*firstneigh;

		ecoul = 0.0;
		if (eflag \|\| vflag) ev_setup(eflag,vflag);
		else evflag = vflag_fdotr = 0;

		double **x = atom->x;
		double **f = atom->f;
		double *q = atom->q;
		int *type = atom->type;
		int nlocal = atom->nlocal;
		int nall = nlocal + atom->nghost;
		double *special_coul = force->special_coul;
		int newton_pair = force->newton_pair;
		double qqrd2e = force->qqrd2e;

		inum = list->inum;
		ilist = list->ilist;
		numneigh = list->numneigh;
		firstneigh = list->firstneigh;

		// loop over neighbors of my atoms

		for (ii = 0; ii < inum; ii++) {
		i = ilist[ii];
		qtmp = q[i];
		xtmp = x[i][0];
		ytmp = x[i][1];
		ztmp = x[i][2];
		itype = type[i];
		jlist = firstneigh[i];
		jnum = numneigh[i];

		for (jj = 0; jj < jnum; jj++) {
		j = jlist[jj];

		if (j < nall) factor_coul = 1.0;
		else {
		factor_coul = special_coul[j/nall];
		j %= nall;
		}

		delx = xtmp - x[j][0];
		dely = ytmp - x[j][1];
		delz = ztmp - x[j][2];
		rsq = delxdelx + delydely + delz*delz;
		jtype = type[j];

		if (rsq < cutsq[itype][jtype]) {
		r2inv = 1.0/rsq;
		r = sqrt(rsq);
		rinv = 1.0/r;
		screening = exp(-kappa*r);
		forcecoul = qqrd2e * qtmpq[j] screening * (kappa + rinv);
		fpair = factor_coulforcecoul r2inv;

		f[i][0] += delx*fpair;
		f[i][1] += dely*fpair;
		f[i][2] += delz*fpair;
		if (newton_pair \|\| j < nlocal) {
		f[j][0] -= delx*fpair;
		f[j][1] -= dely*fpair;
		f[j][2] -= delz*fpair;
		}

		if (eflag) ecoul = factor_coul * qqrd2e * qtmpq[j] rinv * screening;

		if (evflag) ev_tally(i,j,nlocal,newton_pair,
		0.0,ecoul,fpair,delx,dely,delz);
		}
		}
		}

		if (vflag_fdotr) virial_compute();
		}

		/* ----------------------------------------------------------------------
		global settings
		------------------------------------------------------------------------- */

		void PairCoulDebye::settings(int narg, char **arg)
		{
		if (narg != 2) error->all("Illegal pair_style command");

		kappa = atof(arg[0]);
		cut_global = atof(arg[1]);

		// reset cutoffs that have been explicitly set

		if (allocated) {
		int i,j;
		for (i = 1; i <= atom->ntypes; i++)
		for (j = i+1; j <= atom->ntypes; j++)
		if (setflag[i][j]) cut[i][j] = cut_global;
		}
		}

		/* ----------------------------------------------------------------------
		proc 0 writes to restart file
		------------------------------------------------------------------------- */

		void PairCoulDebye::write_restart_settings(FILE *fp)
		{
		fwrite(&cut_global,sizeof(double),1,fp);
		fwrite(&kappa,sizeof(double),1,fp);
		fwrite(&offset_flag,sizeof(int),1,fp);
		fwrite(&mix_flag,sizeof(int),1,fp);
		}

		/* ----------------------------------------------------------------------
		proc 0 reads from restart file, bcasts
		------------------------------------------------------------------------- */

		void PairCoulDebye::read_restart_settings(FILE *fp)
		{
		if (comm->me == 0) {
		fread(&cut_global,sizeof(double),1,fp);
		fread(&kappa,sizeof(double),1,fp);
		fread(&offset_flag,sizeof(int),1,fp);
		fread(&mix_flag,sizeof(int),1,fp);
		}
		MPI_Bcast(&cut_global,1,MPI_DOUBLE,0,world);
		MPI_Bcast(&kappa,1,MPI_DOUBLE,0,world);
		MPI_Bcast(&offset_flag,1,MPI_INT,0,world);
		MPI_Bcast(&mix_flag,1,MPI_INT,0,world);
		}

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

		double PairCoulDebye::single(int i, int j, int itype, int jtype,
		double rsq, double factor_coul, double factor_lj,
		double &fforce)
		{
		double r2inv,r,rinv,forcecoul,phicoul,screening;

		r2inv = 1.0/rsq;
		r = sqrt(rsq);
		rinv = 1.0/r;
		screening = exp(-kappa*r);
		forcecoul = force->qqrd2e * atom->q[i]atom->q[j]
		screening * (kappa + rinv);
		fforce = factor_coulforcecoul r2inv;

		phicoul = force->qqrd2e * atom->q[i]atom->q[j] rinv * screening;
		return factor_coul*phicoul;
		}

src/pair_coul_debye.h

0 → 100644

+37 −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.
		------------------------------------------------------------------------- */

		#ifndef PAIR_COUL_DEBYE_H
		#define PAIR_COUL_DEBYE_H

		#include "pair_coul_cut.h"

		namespace LAMMPS_NS {

		class PairCoulDebye : public PairCoulCut {
		public:
		PairCoulDebye(class LAMMPS *);
		~PairCoulDebye();
		void compute(int, int);
		void settings(int, char **);
		void write_restart_settings(FILE *);
		void read_restart_settings(FILE *);
		double single(int, int, int, int, double, double, double, double &);

		private:
		double kappa;
		};

		}

		#endif

src/pair_lj_cut_coul_cut.h

+2 −2

Original line number	Diff line number	Diff line
		@@ -27,8 +27,8 @@ class PairLJCutCoulCut : public Pair {
		void coeff(int, char **);
		void init_style();
		double init_one(int, int);
		virtual void write_restart(FILE *);
		virtual void read_restart(FILE *);
		void write_restart(FILE *);
		void read_restart(FILE *);
		virtual void write_restart_settings(FILE *);
		virtual void read_restart_settings(FILE *);
		virtual double single(int, int, int, int, double, double, double, double &);

src/pair_lj_cut_coul_debye.cpp

+0 −52

Original line number	Diff line number	Diff line
		@@ -160,58 +160,6 @@ void PairLJCutCoulDebye::settings(int narg, char **arg)
		proc 0 writes to restart file
		------------------------------------------------------------------------- */

		void PairLJCutCoulDebye::write_restart(FILE *fp)
		{
		write_restart_settings(fp);

		int i,j;
		for (i = 1; i <= atom->ntypes; i++)
		for (j = i; j <= atom->ntypes; j++) {
		fwrite(&setflag[i][j],sizeof(int),1,fp);
		if (setflag[i][j]) {
		fwrite(&epsilon[i][j],sizeof(double),1,fp);
		fwrite(&sigma[i][j],sizeof(double),1,fp);
		fwrite(&cut_lj[i][j],sizeof(double),1,fp);
		fwrite(&cut_coul[i][j],sizeof(double),1,fp);
		}
		}
		}

		/* ----------------------------------------------------------------------
		proc 0 reads from restart file, bcasts
		------------------------------------------------------------------------- */

		void PairLJCutCoulDebye::read_restart(FILE *fp)
		{
		read_restart_settings(fp);

		allocate();

		int i,j;
		int me = comm->me;
		for (i = 1; i <= atom->ntypes; i++)
		for (j = i; j <= atom->ntypes; j++) {
		if (me == 0) fread(&setflag[i][j],sizeof(int),1,fp);
		MPI_Bcast(&setflag[i][j],1,MPI_INT,0,world);
		if (setflag[i][j]) {
		if (me == 0) {
		fread(&epsilon[i][j],sizeof(double),1,fp);
		fread(&sigma[i][j],sizeof(double),1,fp);
		fread(&cut_lj[i][j],sizeof(double),1,fp);
		fread(&cut_coul[i][j],sizeof(double),1,fp);
		}
		MPI_Bcast(&epsilon[i][j],1,MPI_DOUBLE,0,world);
		MPI_Bcast(&sigma[i][j],1,MPI_DOUBLE,0,world);
		MPI_Bcast(&cut_lj[i][j],1,MPI_DOUBLE,0,world);
		MPI_Bcast(&cut_coul[i][j],1,MPI_DOUBLE,0,world);
		}
		}
		}

		/* ----------------------------------------------------------------------
		proc 0 writes to restart file
		------------------------------------------------------------------------- */

		void PairLJCutCoulDebye::write_restart_settings(FILE *fp)
		{
		fwrite(&cut_lj_global,sizeof(double),1,fp);

Admin message