From 673202d05f2b0732d2c1be216956a9afa7d11a2f Mon Sep 17 00:00:00 2001
From: Julien Devemy <julien.devemy@uca.fr>
Date: Thu, 14 Jun 2018 13:44:55 +0200
Subject: [PATCH] Bugfix for pair_nm and tail correction
---
src/MISC/pair_nm_cut.cpp | 18 ++++++++++--------
src/MISC/pair_nm_cut_coul_cut.cpp | 17 +++++++++--------
src/MISC/pair_nm_cut_coul_long.cpp | 18 ++++++++++--------
3 files changed, 29 insertions(+), 24 deletions(-)
diff --git a/src/MISC/pair_nm_cut.cpp b/src/MISC/pair_nm_cut.cpp
index 83d734de0f..49ab75194e 100644
--- a/src/MISC/pair_nm_cut.cpp
+++ b/src/MISC/pair_nm_cut.cpp
@@ -15,10 +15,10 @@
Contributing Author: Julien Devemy (ICCF)
------------------------------------------------------------------------- */
-#include <cmath>
-#include <cstdio>
-#include <cstdlib>
-#include <cstring>
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
#include "pair_nm_cut.h"
#include "atom.h"
#include "comm.h"
@@ -282,10 +282,12 @@ double PairNMCut::init_one(int i, int j)
double rrr1 = pow(r0[i][j],nn[i][j])*(1-nn[i][j]);
double rrr2 = pow(r0[i][j],mm[i][j])*(1-mm[i][j]);
- etail_ij = 2.0*MY_PI*all[0]*all[1]*e0nm[i][j] *
- (rr1*pow(cut[i][j],p1)-rr2*pow(cut[i][j],p2));
- ptail_ij = 2.0*MY_PI*all[0]*all[1]*e0nm[i][j] *
- nn[i][j]*mm[i][j]*(rrr1*pow(cut[i][j],p1)-rrr2*pow(cut[i][j],p2));
+ double cut3 = cut[i][j]*cut[i][j]*cut[i][j];
+ ptail_ij = 2.*MY_PI/3.*all[0]*all[1]*e0nm[i][j]*nm[i][j]*cut3 *
+ (pow(r0[i][j]/cut[i][j],nn[i][j])/(nn[i][j]-3) - pow(r0[i][j]/cut[i][j],mm[i][j])/(mm[i][j]-3));
+ etail_ij = 2.*MY_PI*all[0]*all[1]*e0nm[i][j]*cut3 *
+ (mm[i][j]*pow(r0[i][j]/cut[i][j],nn[i][j])/(nn[i][j]-3) - nn[i][j]*pow(r0[i][j]/cut[i][j],mm[i][j])/(mm[i][j]-3));
+
}
return cut[i][j];
diff --git a/src/MISC/pair_nm_cut_coul_cut.cpp b/src/MISC/pair_nm_cut_coul_cut.cpp
index 0c9a43f4dc..cc3a2d8f7d 100644
--- a/src/MISC/pair_nm_cut_coul_cut.cpp
+++ b/src/MISC/pair_nm_cut_coul_cut.cpp
@@ -15,10 +15,10 @@
Contributing Author: Julien Devemy (ICCF)
------------------------------------------------------------------------- */
-#include <cmath>
-#include <cstdio>
-#include <cstdlib>
-#include <cstring>
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
#include "pair_nm_cut_coul_cut.h"
#include "atom.h"
#include "comm.h"
@@ -332,10 +332,11 @@ double PairNMCutCoulCut::init_one(int i, int j)
double rrr1 = pow(r0[i][j],nn[i][j])*(1-nn[i][j]);
double rrr2 = pow(r0[i][j],mm[i][j])*(1-mm[i][j]);
- etail_ij = 2.0*MY_PI*all[0]*all[1]*e0nm[i][j] *
- (rr1*pow(cut_lj[i][j],p1)-rr2*pow(cut_lj[i][j],p2));
- ptail_ij = 2.0*MY_PI*all[0]*all[1]*e0nm[i][j] *
- nn[i][j]*mm[i][j]*(rrr1*pow(cut_lj[i][j],p1)-rrr2*pow(cut_lj[i][j],p2));
+ double cut_lj3 = cut_lj[i][j]*cut_lj[i][j]*cut_lj[i][j];
+ ptail_ij = 2.*MY_PI/3.*all[0]*all[1]*e0nm[i][j]*nm[i][j]*cut_lj3 *
+ (pow(r0[i][j]/cut_lj[i][j],nn[i][j])/(nn[i][j]-3) - pow(r0[i][j]/cut_lj[i][j],mm[i][j])/(mm[i][j]-3));
+ etail_ij = 2.*MY_PI*all[0]*all[1]*e0nm[i][j]*cut_lj3 *
+ (mm[i][j]*pow(r0[i][j]/cut_lj[i][j],nn[i][j])/(nn[i][j]-3) - nn[i][j]*pow(r0[i][j]/cut_lj[i][j],mm[i][j])/(mm[i][j]-3));
}
diff --git a/src/MISC/pair_nm_cut_coul_long.cpp b/src/MISC/pair_nm_cut_coul_long.cpp
index 3296f8fd92..a0aaae5f64 100644
--- a/src/MISC/pair_nm_cut_coul_long.cpp
+++ b/src/MISC/pair_nm_cut_coul_long.cpp
@@ -15,10 +15,10 @@
Contributing Author: Julien Devemy (ICCF)
------------------------------------------------------------------------- */
-#include <cmath>
-#include <cstdio>
-#include <cstdlib>
-#include <cstring>
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
#include "pair_nm_cut_coul_long.h"
#include "atom.h"
#include "comm.h"
@@ -379,10 +379,12 @@ double PairNMCutCoulLong::init_one(int i, int j)
double rrr1 = pow(r0[i][j],nn[i][j])*(1-nn[i][j]);
double rrr2 = pow(r0[i][j],mm[i][j])*(1-mm[i][j]);
- etail_ij = 2.0*MY_PI*all[0]*all[1]*e0nm[i][j] *
- (rr1*pow(cut_lj[i][j],p1)-rr2*pow(cut_lj[i][j],p2));
- ptail_ij = 2.0*MY_PI*all[0]*all[1]*e0nm[i][j]*nn[i][j]*mm[i][j] *
- (rrr1*pow(cut_lj[i][j],p1)-rrr2*pow(cut_lj[i][j],p2));
+ double cut_lj3 = cut_lj[i][j]*cut_lj[i][j]*cut_lj[i][j];
+ ptail_ij = 2.*MY_PI/3.*all[0]*all[1]*e0nm[i][j]*nm[i][j]*cut_lj3 *
+ (pow(r0[i][j]/cut_lj[i][j],nn[i][j])/(nn[i][j]-3) - pow(r0[i][j]/cut_lj[i][j],mm[i][j])/(mm[i][j]-3));
+ etail_ij = 2.*MY_PI*all[0]*all[1]*e0nm[i][j]*cut_lj3 *
+ (mm[i][j]*pow(r0[i][j]/cut_lj[i][j],nn[i][j])/(nn[i][j]-3) - nn[i][j]*pow(r0[i][j]/cut_lj[i][j],mm[i][j])/(mm[i][j]-3));
+
}
return cut;
--
GitLab