diff --git a/src/MANYBODY/pair_eam.cpp b/src/MANYBODY/pair_eam.cpp
index 3f6972bc5e..5e5f8f3930 100644
--- a/src/MANYBODY/pair_eam.cpp
+++ b/src/MANYBODY/pair_eam.cpp
@@ -206,6 +206,8 @@ void PairEAM::compute(int eflag, int vflag)
 
   // fp = derivative of embedding energy at each atom
   // phi = embedding energy at each atom
+  // if rho > rhomax (e.g. due to close approach of two atoms),
+  //   will exceed table, so add linear term to conserve energy
 
   for (ii = 0; ii < inum; ii++) {
     i = ilist[ii];
@@ -218,6 +220,7 @@ void PairEAM::compute(int eflag, int vflag)
     fp[i] = (coeff[0]*p + coeff[1])*p + coeff[2];
     if (eflag) {
       phi = ((coeff[3]*p + coeff[4])*p + coeff[5])*p + coeff[6];
+      if (rho[i] > rhomax) phi += fp[i] * (rho[i]-rhomax);
       if (eflag_global) eng_vdwl += phi;
       if (eflag_atom) eatom[i] += phi;
     }
@@ -491,7 +494,7 @@ void PairEAM::file2array()
   // active means some element is pointing at it via map
 
   int active;
-  double rmax,rhomax;
+  double rmax;
   dr = drho = rmax = rhomax = 0.0;
 
   for (int i = 0; i < nfuncfl; i++) {
diff --git a/src/MANYBODY/pair_eam.h b/src/MANYBODY/pair_eam.h
index e91e73b5d5..96c36b99d2 100644
--- a/src/MANYBODY/pair_eam.h
+++ b/src/MANYBODY/pair_eam.h
@@ -41,7 +41,7 @@ class PairEAM : public Pair {
 
   // potentials in spline form used for force computation
 
-  double dr,rdr,drho,rdrho;
+  double dr,rdr,drho,rdrho,rhomax;
   double ***rhor_spline,***frho_spline,***z2r_spline;
 
   PairEAM(class LAMMPS *);
diff --git a/src/MANYBODY/pair_eam_alloy.cpp b/src/MANYBODY/pair_eam_alloy.cpp
index 0f926a8cb7..318d6421e6 100644
--- a/src/MANYBODY/pair_eam_alloy.cpp
+++ b/src/MANYBODY/pair_eam_alloy.cpp
@@ -219,6 +219,7 @@ void PairEAMAlloy::file2array()
   nr = setfl->nr;
   drho = setfl->drho;
   dr = setfl->dr;
+  rhomax = (nrho-1) * drho;
 
   // ------------------------------------------------------------------
   // setup frho arrays
diff --git a/src/MANYBODY/pair_eam_fs.cpp b/src/MANYBODY/pair_eam_fs.cpp
index 29f787b296..875822fd5c 100644
--- a/src/MANYBODY/pair_eam_fs.cpp
+++ b/src/MANYBODY/pair_eam_fs.cpp
@@ -224,6 +224,7 @@ void PairEAMFS::file2array()
   nr = fs->nr;
   drho = fs->drho;
   dr = fs->dr;
+  rhomax = (nrho-1) * drho;
 
   // ------------------------------------------------------------------
   // setup frho arrays
diff --git a/src/OPT/pair_eam_opt.cpp b/src/OPT/pair_eam_opt.cpp
index c0cb4165ec..ba0512df39 100644
--- a/src/OPT/pair_eam_opt.cpp
+++ b/src/OPT/pair_eam_opt.cpp
@@ -233,6 +233,8 @@ void PairEAMOpt::eval()
 
   // fp = derivative of embedding energy at each atom
   // phi = embedding energy at each atom
+  // if rho > rhomax (e.g. due to close approach of two atoms),
+  //   will exceed table, so add linear term to conserve energy
 
   for (ii = 0; ii < inum; ii++) {
     i = ilist[ii];
@@ -244,6 +246,7 @@ void PairEAMOpt::eval()
     fp[i] = (coeff[0]*p + coeff[1])*p + coeff[2];
     if (EFLAG) {
       double phi = ((coeff[3]*p + coeff[4])*p + coeff[5])*p + coeff[6];
+      if (rho[i] > rhomax) phi += fp[i] * (rho[i]-rhomax);
       if (eflag_global) eng_vdwl += phi;
       if (eflag_atom) eatom[i] += phi;
     }