From 405c6a0e2ca0e7fef785a35d47786af8155fba66 Mon Sep 17 00:00:00 2001
From: Steve Plimpton <sjplimp@gmail.com>
Date: Mon, 1 May 2023 14:41:23 -0600
Subject: [PATCH 1/2] code and documentation

---
 doc/fix_emit_surf.html |  40 ++++++++++++----
 doc/fix_emit_surf.txt  |  40 ++++++++++++----
 src/fix_emit_surf.cpp  | 105 ++++++++++++++++++++++++++++-------------
 src/fix_emit_surf.h    |   6 ++-
 4 files changed, 136 insertions(+), 55 deletions(-)

diff --git a/doc/fix_emit_surf.html b/doc/fix_emit_surf.html
index 05d9bed0c..275a150ce 100644
--- a/doc/fix_emit_surf.html
+++ b/doc/fix_emit_surf.html
@@ -28,6 +28,7 @@ <H3>fix emit/surf command
 <LI>keyword = <I>n</I> or <I>normal</I> or <I>nevery</I> or <I>perspecies</I> or <I>region</I> or <I>subsonic</I> 
 
 <PRE>  <I>n</I> value = Np = number of particles to create
+                   Np can be a variable (see below)
   <I>normal</I> value = yes or no = emit normal to surface elements or with streaming velocity
   <I>nevery</I> value = Nstep = add particles every this many timesteps
   <I>perspecies</I> value = <I>yes</I> or <I>no</I>
@@ -71,9 +72,11 @@ <H3>fix emit/surf command
 within a grid cell is based on this flux and additional global, flow,
 and surface element properties:
 </P>
-<UL><LI>global property: <I>fnum</I> ratio as specified by the <A HREF = "global.html"">global</A> command
-<LI>flow properties: number density, streaming velocity, and thermal temperature
-<LI>surface element properties: portion of surface element area that overlaps with the grid cell and its orientation relative to the streaming velocity 
+<UL><LI>global property: <I>fnum</I> ratio as specified by the
+<LI><A HREF = "global.html"">global</A> command flow properties: number density,
+<LI>streaming velocity, and thermal temperature surface element
+<LI>properties: portion of surface element area that overlaps with the
+<LI>grid cell and its orientation relative to the streaming velocity 
 </UL>
 <P>The flow properties are defined for the specified mixture via the
 <A HREF = "mixture.html">mixture</A> command.
@@ -116,11 +119,28 @@ <H3>fix emit/surf command
 <P>The <I>n</I> keyword can alter how many particles are added, which can be
 useful for debugging purposes.  If <I>Np</I> is set to 0, then the number
 of added particles is a function of <I>fnum</I>, <I>nrho</I>, and other mixture
-settings, as described above.  If <I>Np</I> is set to a value > 0, then the
-<I>fnum</I> and <I>nrho</I> settings are ignored, and exactly <I>Np</I> particles are
-added on each insertion timestep.  This is done by dividing <I>Np</I> by
-the total number of grid cell/surface element pairs and adding an
-equal number of particles per pair.
+settings, as described above.
+</P>
+<P>If <I>Np</I> is set to a value > 0, then the <I>fnum</I> and <I>nrho</I> settings are
+ignored, and roughly <I>Np</I> particles are added on each insertion
+timestep.  For each grid cell/surface element pair, its target number
+of emitted particles is set to its fraction of the total emission area
+(for all grid cell/surface element pairs), multiplied by <I>Np</I>.  If
+that results in a fractional value, then an extra particle is emitted
+depending on the value of a random number, as explained above.
+</P>
+<P>The <I>Np</I> value can be also be specified as an equal-style
+<A HREF = "variable.html">variable</A>.  If the value is a variable, it should be
+specified as v_name, where name is the variable name.  In this case,
+the variable will be evaluated on each emission timestep, and its
+value used as <I>Np</I> on that step to determine the target number of
+emitted particles for each grid cell/surface element pair, the same as
+described in the preceeding paragraph.
+</P>
+<P>Equal-style variables can specify formulas with various mathematical
+functions, and include <A HREF = "status_style.html">stats_style</A> command
+keywords for the simulation box parameters and timestep and elapsed
+time.  Thus it is easy to specify a time-dependent value of <I>Np</I>.
 </P>
 <P>The <I>normal</I> keyword can be used to alter how velocities are set for
 added particles.  If <I>normal</I> is set to <I>no</I>, then a particle's
@@ -240,8 +260,8 @@ <H3>fix emit/surf command
 </P>
 <P><B>Restrictions:</B>
 </P>
-<P>A <I>n</I> setting of <I>Np</I> > 0 can only be used with a <I>perspecies</I> setting
-of <I>no</I>.
+<P>A <I>n</I> setting of <I>Np</I> > 0 or <I>Np</I> as a variable can only be used with
+a <I>perspecies</I> setting of <I>no</I>.
 </P>
 <P>If <I>normal</I> is set to <I>no</I>, which is the default, then unlike the <A HREF = "fix_emit/face.html">fix
 emit/face</A> command, no warning is issued if a
diff --git a/doc/fix_emit_surf.txt b/doc/fix_emit_surf.txt
index e27f5aed2..dee5e80e5 100644
--- a/doc/fix_emit_surf.txt
+++ b/doc/fix_emit_surf.txt
@@ -19,6 +19,7 @@ group-ID = ID of surface group that emits particles :l
 zero or more keyword/value pairs may be appended :l
 keyword = {n} or {normal} or {nevery} or {perspecies} or {region} or {subsonic} :l
   {n} value = Np = number of particles to create
+                   Np can be a variable (see below)
   {normal} value = yes or no = emit normal to surface elements or with streaming velocity
   {nevery} value = Nstep = add particles every this many timesteps
   {perspecies} value = {yes} or {no}
@@ -62,9 +63,11 @@ given by equation 4.22 of "(Bird94)"_#Bird94.  The number of particles
 within a grid cell is based on this flux and additional global, flow,
 and surface element properties:
 
-global property: {fnum} ratio as specified by the "global"_global.html" command
-flow properties: number density, streaming velocity, and thermal temperature
-surface element properties: portion of surface element area that overlaps with the grid cell and its orientation relative to the streaming velocity :ul
+global property: {fnum} ratio as specified by the
+"global"_global.html" command flow properties: number density,
+streaming velocity, and thermal temperature surface element
+properties: portion of surface element area that overlaps with the
+grid cell and its orientation relative to the streaming velocity :ul
 
 The flow properties are defined for the specified mixture via the
 "mixture"_mixture.html command.
@@ -107,11 +110,28 @@ of particles emitting from each surface element.
 The {n} keyword can alter how many particles are added, which can be
 useful for debugging purposes.  If {Np} is set to 0, then the number
 of added particles is a function of {fnum}, {nrho}, and other mixture
-settings, as described above.  If {Np} is set to a value > 0, then the
-{fnum} and {nrho} settings are ignored, and exactly {Np} particles are
-added on each insertion timestep.  This is done by dividing {Np} by
-the total number of grid cell/surface element pairs and adding an
-equal number of particles per pair.
+settings, as described above.
+
+If {Np} is set to a value > 0, then the {fnum} and {nrho} settings are
+ignored, and roughly {Np} particles are added on each insertion
+timestep.  For each grid cell/surface element pair, its target number
+of emitted particles is set to its fraction of the total emission area
+(for all grid cell/surface element pairs), multiplied by {Np}.  If
+that results in a fractional value, then an extra particle is emitted
+depending on the value of a random number, as explained above.
+
+The {Np} value can be also be specified as an equal-style
+"variable"_variable.html.  If the value is a variable, it should be
+specified as v_name, where name is the variable name.  In this case,
+the variable will be evaluated on each emission timestep, and its
+value used as {Np} on that step to determine the target number of
+emitted particles for each grid cell/surface element pair, the same as
+described in the preceeding paragraph.
+
+Equal-style variables can specify formulas with various mathematical
+functions, and include "stats_style"_status_style.html command
+keywords for the simulation box parameters and timestep and elapsed
+time.  Thus it is easy to specify a time-dependent value of {Np}.
 
 The {normal} keyword can be used to alter how velocities are set for
 added particles.  If {normal} is set to {no}, then a particle's
@@ -231,8 +251,8 @@ a run is performed.
 
 [Restrictions:]
 
-A {n} setting of {Np} > 0 can only be used with a {perspecies} setting
-of {no}.
+A {n} setting of {Np} > 0 or {Np} as a variable can only be used with
+a {perspecies} setting of {no}.
 
 If {normal} is set to {no}, which is the default, then unlike the "fix
 emit/face"_fix_emit/face.html command, no warning is issued if a
diff --git a/src/fix_emit_surf.cpp b/src/fix_emit_surf.cpp
index 3201608e4..1e70f35eb 100644
--- a/src/fix_emit_surf.cpp
+++ b/src/fix_emit_surf.cpp
@@ -25,6 +25,7 @@
 #include "cut2d.h"
 #include "cut3d.h"
 #include "input.h"
+#include "variable.h"
 #include "comm.h"
 #include "random_knuth.h"
 #include "math_extra.h"
@@ -37,6 +38,7 @@ using namespace MathConst;
 
 enum{PKEEP,PINSERT,PDONE,PDISCARD,PENTRY,PEXIT,PSURF};   // several files
 enum{NOSUBSONIC,PTBOTH,PONLY};
+enum{FLOW,CONSTANT,VARIABLE};
 
 #define DELTATASK 256
 #define TEMPLIMIT 1.0e5
@@ -60,6 +62,8 @@ FixEmitSurf::FixEmitSurf(SPARTA *sparta, int narg, char **arg) :
   // optional args
 
   np = 0;
+  npmode = FLOW;
+  npstr = NULL;
   normalflag = 0;
   subsonic = 0;
   subsonic_style = NOSUBSONIC;
@@ -73,9 +77,10 @@ FixEmitSurf::FixEmitSurf(SPARTA *sparta, int narg, char **arg) :
   if (!surf->exist)
     error->all(FLERR,"Fix emit/surf requires surface elements");
   if (surf->implicit)
-    error->all(FLERR,"Fix emit/surf not allowed for implicits surfaces");
-  if (np > 0 && perspecies)
-    error->all(FLERR,"Cannot use fix emit/face n > 0 with perspecies yes");
+    error->all(FLERR,"Fix emit/surf not allowed for implicit surfaces");
+  if ((npmode == CONSTANT || npmode == VARIABLE) && perspecies)
+    error->all(FLERR,"Cannot use fix emit/surf with n a constant or variable "
+               "with perspecies yes");
 
   // task list and subsonic data structs
 
@@ -97,6 +102,8 @@ FixEmitSurf::FixEmitSurf(SPARTA *sparta, int narg, char **arg) :
 
 FixEmitSurf::~FixEmitSurf()
 {
+  delete [] npstr;
+  
   for (int i = 0; i < ntaskmax; i++) {
     delete [] tasks[i].ntargetsp;
     delete [] tasks[i].vscale;
@@ -177,6 +184,16 @@ void FixEmitSurf::init()
     }
   }
 
+  // check variable for npmode = VARIABLE
+
+  if (npmode == VARIABLE) {
+    npvar = input->variable->find(npstr);
+    if (npvar < 0)
+      error->all(FLERR,"Fix emit/surf variable name does not exist");
+    if (!input->variable->equal_style(npvar))
+      error->all(FLERR,"Fix emit/surf  variable is invalid style");
+  }
+  
   // create tasks for all grid cells
 
   grid_changed();
@@ -193,24 +210,19 @@ void FixEmitSurf::grid_changed()
 {
   create_tasks();
 
-  // if Np > 0, nper = # of insertions per task
-  // set nthresh so as to achieve exactly Np insertions
-  // tasks > tasks_with_no_extra need to insert 1 extra particle
-  // NOTE: setting same # of insertions per task
-  //       should weight by overlap area of cell/surf
-
-  if (np > 0) {
-    int all,nupto,tasks_with_no_extra;
-    MPI_Allreduce(&ntask,&all,1,MPI_INT,MPI_SUM,world);
-    if (all) {
-      npertask = np / all;
-      tasks_with_no_extra = all - (np % all);
-    } else npertask = tasks_with_no_extra = 0;
-
-    MPI_Scan(&ntask,&nupto,1,MPI_INT,MPI_SUM,world);
-    if (tasks_with_no_extra < nupto-ntask) nthresh = 0;
-    else if (tasks_with_no_extra >= nupto) nthresh = ntask;
-    else nthresh = tasks_with_no_extra - (nupto-ntask);
+  // for MODE = CONSTANT or VARIABLE
+  // set per-task ntarget to fraction of its area / total area
+
+  if (npmode != FLOW) {
+    double areasum_me = 0.0;
+    for (int i = 0; i < ntask; i++)
+      areasum_me += tasks[i].area;
+
+    double areasum;
+    MPI_Allreduce(&areasum_me,&areasum,1,MPI_DOUBLE,MPI_SUM,world);
+
+    for (int i = 0; i < ntask; i++)
+      tasks[i].ntarget = tasks[i].area / areasum;
   }
 }
 
@@ -219,9 +231,6 @@ void FixEmitSurf::grid_changed()
 void FixEmitSurf::setup()
 {
   // needed for Kokkos because pointers are changed in UpdateKokkos::setup()
-
-  //lines = surf->lines;
-  //tris = surf->tris;
 }
 
 /* ----------------------------------------------------------------------
@@ -383,6 +392,7 @@ void FixEmitSurf::create_task(int icell)
     }
 
     // set ntarget and ntargetsp via mol_inflow()
+    // will be overwritten if mode != FLOW
     // skip task if final ntarget = 0.0, due to large outbound vstream
     // do not skip for subsonic since it resets ntarget every step
 
@@ -441,6 +451,14 @@ void FixEmitSurf::perform_task()
 
   if (subsonic) subsonic_inflow();
 
+  // if npmode = VARIABLE, set npcurrent to variable evaluation
+
+  double npcurrent;
+  if (npmode == VARIABLE) {
+    npcurrent = input->variable->compute_equal(npvar);
+    if (npcurrent <= 0.0) error->all(FLERR,"Fix emit/surf Np <= 0.0");
+  }
+  
   // insert particles for each task = cell/surf pair
   // ntarget/ninsert is either perspecies or for all species
   // for one particle:
@@ -449,7 +467,7 @@ void FixEmitSurf::perform_task()
   //       if normalflag, mag of vstream is applied to surf normal dir
   //       first stage: normal dimension (normal)
   //       second stage: parallel dimensions (tan1,tan2)
-
+  
   // double while loop until randomized particle velocity meets 2 criteria
   // inner do-while loop:
   //   v = vstream-component + vthermal is into simulation box
@@ -458,6 +476,7 @@ void FixEmitSurf::perform_task()
   //   shift Maxwellian distribution by stream velocity component
   //   see Bird 1994, p 259, eq 12.5
 
+
   Surf::Line *lines = surf->lines;
   Surf::Tri *tris = surf->tris;
 
@@ -483,6 +502,8 @@ void FixEmitSurf::perform_task()
     if (!normalflag) indot = vstream[0]*normal[0] + vstream[1]*normal[1] +
                        vstream[2]*normal[2];
 
+    // perspecies yes
+    
     if (perspecies) {
       for (isp = 0; isp < nspecies; isp++) {
         ispecies = species[isp];
@@ -567,14 +588,20 @@ void FixEmitSurf::perform_task()
         nsingle += nactual;
       }
 
+    // perspecies no
+
     } else {
-      if (np == 0) {
-        ntarget = tasks[i].ntarget+random->uniform();
-        ninsert = static_cast<int> (ntarget);
-      } else {
-        ninsert = npertask;
-        if (i >= nthresh) ninsert++;
-      }
+
+      // set ntarget for insertion mode FLOW, CONSTANT, or VARIABLE
+      // for FLOW: ntarget is already set within task
+      // for CONSTANT or VARIABLE: task narget is fraction of its surf's area
+      //   scale fraction by np or npcurrent (variable evaluation)
+      // ninsert = rounded-down (ntarget + random number)
+      
+      if (npmode == FLOW) ntarget = tasks[i].ntarget;
+      else if (npmode == CONSTANT) ntarget = np * tasks[i].ntarget;
+      else if (npmode == VARIABLE) ntarget = npcurrent * tasks[i].ntarget;
+      ninsert = static_cast<int> (ntarget + random->uniform());
 
       nactual = 0;
       for (int m = 0; m < ninsert; m++) {
@@ -965,8 +992,18 @@ int FixEmitSurf::option(int narg, char **arg)
 {
   if (strcmp(arg[0],"n") == 0) {
     if (2 > narg) error->all(FLERR,"Illegal fix emit/]surf/normal command");
-    np = atoi(arg[1]);
-    if (np <= 0) error->all(FLERR,"Illegal fix emit/surf/normal command");
+
+    if (strstr(arg[1],"v_") == arg[1]) {
+      npmode = VARIABLE;
+      int n = strlen(&arg[1][2]) + 1;
+      npstr = new char[n];
+      strcpy(npstr,&arg[1][2]);
+    } else {
+      np = atoi(arg[1]);
+      if (np <= 0) error->all(FLERR,"Illegal fix emit/surf/normal command");
+      if (np == 0) npmode = FLOW;
+      else npmode = CONSTANT;
+    }
     return 2;
   }
 
diff --git a/src/fix_emit_surf.h b/src/fix_emit_surf.h
index d7db6f3b6..ebcaff322 100644
--- a/src/fix_emit_surf.h
+++ b/src/fix_emit_surf.h
@@ -37,11 +37,15 @@ class FixEmitSurf : public FixEmit {
   void grid_changed();
 
  private:
-  int imix,groupbit,np,normalflag,subsonic,subsonic_style,subsonic_warning;
+  int imix,groupbit,normalflag,subsonic,subsonic_style,subsonic_warning;
   int npertask,nthresh;
   double psubsonic,tsubsonic,nsubsonic;
   double tprefactor,soundspeed_mixture;
 
+  int npmode,np;    // npmode = FLOW,CONSTANT,VARIABLE
+  int npvar;
+  char *npstr;
+  
   // copies of data from other classes
 
   int dimension,nspecies;

From 58df4df9766898d149d46a18c402321040ee0dab Mon Sep 17 00:00:00 2001
From: Stan Gerald Moore <stamoor@sandia.gov>
Date: Fri, 2 Jun 2023 12:32:40 -0600
Subject: [PATCH 2/2] Small cleanup

---
 src/create_particles.cpp | 60 ++++++++++++++++++++--------------------
 src/cut2d.h              |  6 ++--
 src/cut3d.h              |  6 ++--
 src/fix_ablate.cpp       | 10 +++----
 src/fix_emit_surf.cpp    | 19 ++++---------
 src/fix_emit_surf.h      |  3 +-
 src/fix_grid_check.cpp   |  8 +++---
 src/grid.cpp             |  2 +-
 src/grid_surf.cpp        | 14 +++++-----
 9 files changed, 60 insertions(+), 68 deletions(-)

diff --git a/src/create_particles.cpp b/src/create_particles.cpp
index fd6a66458..2a05c4317 100644
--- a/src/create_particles.cpp
+++ b/src/create_particles.cpp
@@ -319,7 +319,7 @@ void CreateParticles::command(int narg, char **arg)
   double time1 = MPI_Wtime();
 
   bigint nprevious = particle->nglobal;
-  
+
   if (single) create_single();
   else if (!globalflag) {
     if (twopass) create_local_twopass();
@@ -459,22 +459,22 @@ void CreateParticles::create_local()
   // insertvol = subset of flowvol for cells eligible for insertion
   //   insertvol = flowvol if cutflag = 1
   //   insertvol < flowvol possible if cutflag = 0 (no cut cells)
-  
+
   double flowvolme = 0.0;
   double insertvolme = 0.0;
-  
+
   for (int icell = 0; icell < nglocal; icell++) {
     if (cinfo[icell].type == INSIDE) continue;
     if (cells[icell].nsplit > 1) continue;
     if (region && region->bboxflag &&
         outside_region(dimension,cells[icell].lo,cells[icell].hi))
       continue;
-    
+
     flowvolme += cinfo[icell].volume / cinfo[icell].weight;
     if (!cutflag && cells[icell].nsurf) continue;
     insertvolme += cinfo[icell].volume / cinfo[icell].weight;
   }
-  
+
   // calculate total Np if not set explicitly
   // based on total flowvol and mixture density
 
@@ -483,9 +483,9 @@ void CreateParticles::create_local()
     MPI_Allreduce(&flowvolme,&flowvol,1,MPI_DOUBLE,MPI_SUM,world);
     np = particle->mixture[imix]->nrho * flowvol / update->fnum;
   }
-  
+
   // gather cummulative insertion volumes across all procs
-  
+
   double volupto;
   MPI_Scan(&insertvolme,&volupto,1,MPI_DOUBLE,MPI_SUM,world);
 
@@ -542,7 +542,7 @@ void CreateParticles::create_local()
   double x[3],v[3],xcell[3],vstream_variable[3];
   double ntarget,scale,rn,vn,vr,theta1,theta2,erot,evib;
   double *lo,*hi;
-  
+
   double tempscale = 1.0;
   double sqrttempscale = 1.0;
 
@@ -557,9 +557,9 @@ void CreateParticles::create_local()
         outside_region(dimension,cells[icell].lo,cells[icell].hi))
       continue;
     if (!cutflag && cells[icell].nsurf) continue;
-    
+
     volsum += cinfo[icell].volume / cinfo[icell].weight;
-    
+
     ntarget = nme * volsum/insertvolme - nprev;
     npercell = static_cast<int> (ntarget);
 
@@ -568,7 +568,7 @@ void CreateParticles::create_local()
 
     lo = cells[icell].lo;
     hi = cells[icell].hi;
-    
+
     if (densflag) {
       scale = density_variable(lo,hi);
       ntarget *= scale;
@@ -580,7 +580,7 @@ void CreateParticles::create_local()
 
     if (cells[icell].nsurf)
       pflag = grid->point_outside_surfs(icell,xcell);
-    
+
     for (int m = 0; m < ncreate; m++) {
 
       // generate random position X for new particle
@@ -593,7 +593,7 @@ void CreateParticles::create_local()
       // if surfs, check if random position is in flow region
       // if subcell, also check if in correct subcell
       // if not, attempt new insertion up to MAXATTEMPT times
-      
+
       if (cells[icell].nsurf && pflag) {
         int nattempt = 0;
         while (nattempt < MAXATTEMPT) {
@@ -606,7 +606,7 @@ void CreateParticles::create_local()
           }
 
           nattempt++;
-          
+
           x[0] = lo[0] + random->uniform() * (hi[0]-lo[0]);
           x[1] = lo[1] + random->uniform() * (hi[1]-lo[1]);
           x[2] = lo[2] + random->uniform() * (hi[2]-lo[2]);
@@ -614,7 +614,7 @@ void CreateParticles::create_local()
         }
 
         // particle insertion was unsuccessful
-        
+
         if (nattempt >= MAXATTEMPT) continue;
       }
 
@@ -662,7 +662,7 @@ void CreateParticles::create_local()
       evib = particle->evib(ispecies,temp_vib*tempscale,random);
 
       id = MAXSMALLINT*random->uniform();
-      
+
       particle->add_particle(id,ispecies,icell,x,v,erot,evib);
 
       if (nfix_update_custom)
@@ -708,17 +708,17 @@ void CreateParticles::create_local_twopass()
   // insertvol = subset of flowvol for cells eligible for insertion
   //   insertvol = flowvol if cutflag = 1
   //   insertvol < flowvol possible if cutflag = 0 (no cut cells)
-  
+
   double flowvolme = 0.0;
   double insertvolme = 0.0;
-  
+
   for (int icell = 0; icell < nglocal; icell++) {
     if (cinfo[icell].type == INSIDE) continue;
     if (cells[icell].nsplit > 1) continue;
     if (region && region->bboxflag &&
         outside_region(dimension,cells[icell].lo,cells[icell].hi))
       continue;
-    
+
     flowvolme += cinfo[icell].volume / cinfo[icell].weight;
     if (!cutflag && cells[icell].nsurf) continue;
     insertvolme += cinfo[icell].volume / cinfo[icell].weight;
@@ -732,9 +732,9 @@ void CreateParticles::create_local_twopass()
     MPI_Allreduce(&flowvolme,&flowvol,1,MPI_DOUBLE,MPI_SUM,world);
     np = particle->mixture[imix]->nrho * flowvol / update->fnum;
   }
-  
+
   // gather cummulative insertion volumes across all procs
-  
+
   double volupto;
   MPI_Scan(&insertvolme,&volupto,1,MPI_DOUBLE,MPI_SUM,world);
 
@@ -799,7 +799,7 @@ void CreateParticles::create_local_twopass()
 
   // first pass, just calculate # of particles to create
   // ncreate_values[icell] = # of particles to create in ICELL
-  
+
   int *ncreate_values;
   memory->create(ncreate_values, nglocal, "create_particles:ncreate");
 
@@ -813,7 +813,7 @@ void CreateParticles::create_local_twopass()
     if (!cutflag && cells[icell].nsurf) continue;
 
     volsum += cinfo[icell].volume / cinfo[icell].weight;
-    
+
     ntarget = nme * volsum/insertvolme - nprev;
     npercell = static_cast<int> (ntarget);
 
@@ -839,7 +839,7 @@ void CreateParticles::create_local_twopass()
   }
 
   // second pass, create particles using ncreate_values
-  
+
   for (int icell = 0; icell < nglocal; icell++) {
     if (cinfo[icell].type == INSIDE) continue;
     if (cells[icell].nsplit > 1) continue;
@@ -858,7 +858,7 @@ void CreateParticles::create_local_twopass()
 
     if (cells[icell].nsurf)
       pflag = grid->point_outside_surfs(icell,xcell);
-    
+
     for (int m = 0; m < ncreate; m++) {
 
       // generate random position X for new particle
@@ -871,7 +871,7 @@ void CreateParticles::create_local_twopass()
       // if surfs, check if random position is in flow region
       // if subcell, also check if in correct subcell
       // if not, attempt new insertion up to MAXATTEMPT times
-      
+
       if (cells[icell].nsurf && pflag) {
         int nattempt = 0;
         while (nattempt < MAXATTEMPT) {
@@ -884,7 +884,7 @@ void CreateParticles::create_local_twopass()
           }
 
           nattempt++;
-          
+
           x[0] = lo[0] + random->uniform() * (hi[0]-lo[0]);
           x[1] = lo[1] + random->uniform() * (hi[1]-lo[1]);
           x[2] = lo[2] + random->uniform() * (hi[2]-lo[2]);
@@ -892,7 +892,7 @@ void CreateParticles::create_local_twopass()
         }
 
         // particle insertion was unsuccessful
-        
+
         if (nattempt >= MAXATTEMPT) continue;
       }
 
@@ -940,7 +940,7 @@ void CreateParticles::create_local_twopass()
       evib = particle->evib(ispecies,temp_vib*tempscale,random);
 
       id = MAXSMALLINT*random->uniform();
-      
+
       particle->add_particle(id,ispecies,icell,x,v,erot,evib);
 
       if (nfix_update_custom)
@@ -950,7 +950,7 @@ void CreateParticles::create_local_twopass()
   }
 
   // clean up
-  
+
   memory->destroy(ncreate_values);
 
   delete random;
diff --git a/src/cut2d.h b/src/cut2d.h
index 6ee1e81f4..90f03ed03 100644
--- a/src/cut2d.h
+++ b/src/cut2d.h
@@ -81,11 +81,11 @@ class Cut2d : protected Pointers {
   int split(cellint, double *, double *, int, surfint *,
             double *&, int *, int *, int &, double *);
   int split_face(int, int, double *, double *);
-  
+
   int clip_external(double *, double *, double *, double *, double *);
   int sameedge(double *, double *);
   int sameedge_external(double *, double *, double *, double *);
-  
+
  private:
   int axisymmetric;
   int implicit;
@@ -103,7 +103,7 @@ class Cut2d : protected Pointers {
   MyVec<int> used;       // 0/1 flag for each point when walking loops
 
   // methods
-  
+
   void build_clines();
   int weiler_build();
   void weiler_loops();
diff --git a/src/cut3d.h b/src/cut3d.h
index 41721ec28..38701d50b 100644
--- a/src/cut3d.h
+++ b/src/cut3d.h
@@ -31,7 +31,7 @@ class Cut3d : protected Pointers {
   int surf2grid_list(cellint, double *, double *, int, surfint *,
                      surfint *, int);
   int surf2grid_one(double *, double *, double *, double *, double *);
-  
+
   int split(cellint, double *, double *, int, surfint *,
             double *&, int *, int *, int &, double *);
 
@@ -39,7 +39,7 @@ class Cut3d : protected Pointers {
                     double *, double *, double *);
   int sameface(double *, double *, double *);
   int sameface_external(double *, double *, double *, double *, double *);
-  
+
  private:
   int implicit;
 
@@ -120,7 +120,7 @@ class Cut3d : protected Pointers {
   class Cut2d *cut2d;
 
   // methods
-  
+
   int clip(double *, double *, double *);
   int split_try(cellint, int, surfint *,
                 double *&, int *, int *, int &, double *, int &);
diff --git a/src/fix_ablate.cpp b/src/fix_ablate.cpp
index e5592be6b..e56262c66 100644
--- a/src/fix_ablate.cpp
+++ b/src/fix_ablate.cpp
@@ -608,18 +608,18 @@ void FixAblate::create_surfs(int outflag)
     particles[i].flag = PKEEP;
     icell = particles[i].icell;
     if (cells[icell].nsurf == 0) continue;
-    
+
     x = particles[i].x;
 
     // check that particle is outside surfs
     // if no xcell found, cannot check
-    
+
     pflag = grid->point_outside_surfs(icell,xcell);
     if (!pflag) continue;
     pflag = grid->outside_surfs(icell,x,xcell);
-    
+
     // check that particle is in correct split subcell
-    
+
     if (pflag && cells[icell].nsplit <= 0) {
       splitcell = sinfo[cells[icell].isplit].icell;
       if (dim == 2) subcell = update->split2d(splitcell,x);
@@ -628,7 +628,7 @@ void FixAblate::create_surfs(int outflag)
     }
 
     // discard the particle if either test failed
-    
+
     if (!pflag) {
       particles[i].flag = PDISCARD;
       // DEBUG - print message about MC flags for cell of deleted particle
diff --git a/src/fix_emit_surf.cpp b/src/fix_emit_surf.cpp
index 1e70f35eb..e01ec6964 100644
--- a/src/fix_emit_surf.cpp
+++ b/src/fix_emit_surf.cpp
@@ -103,7 +103,7 @@ FixEmitSurf::FixEmitSurf(SPARTA *sparta, int narg, char **arg) :
 FixEmitSurf::~FixEmitSurf()
 {
   delete [] npstr;
-  
+
   for (int i = 0; i < ntaskmax; i++) {
     delete [] tasks[i].ntargetsp;
     delete [] tasks[i].vscale;
@@ -193,7 +193,7 @@ void FixEmitSurf::init()
     if (!input->variable->equal_style(npvar))
       error->all(FLERR,"Fix emit/surf  variable is invalid style");
   }
-  
+
   // create tasks for all grid cells
 
   grid_changed();
@@ -226,13 +226,6 @@ void FixEmitSurf::grid_changed()
   }
 }
 
-/* ---------------------------------------------------------------------- */
-
-void FixEmitSurf::setup()
-{
-  // needed for Kokkos because pointers are changed in UpdateKokkos::setup()
-}
-
 /* ----------------------------------------------------------------------
    create task for one grid cell
    add them to tasks list and increment ntasks
@@ -458,7 +451,7 @@ void FixEmitSurf::perform_task()
     npcurrent = input->variable->compute_equal(npvar);
     if (npcurrent <= 0.0) error->all(FLERR,"Fix emit/surf Np <= 0.0");
   }
-  
+
   // insert particles for each task = cell/surf pair
   // ntarget/ninsert is either perspecies or for all species
   // for one particle:
@@ -467,7 +460,7 @@ void FixEmitSurf::perform_task()
   //       if normalflag, mag of vstream is applied to surf normal dir
   //       first stage: normal dimension (normal)
   //       second stage: parallel dimensions (tan1,tan2)
-  
+
   // double while loop until randomized particle velocity meets 2 criteria
   // inner do-while loop:
   //   v = vstream-component + vthermal is into simulation box
@@ -503,7 +496,7 @@ void FixEmitSurf::perform_task()
                        vstream[2]*normal[2];
 
     // perspecies yes
-    
+
     if (perspecies) {
       for (isp = 0; isp < nspecies; isp++) {
         ispecies = species[isp];
@@ -597,7 +590,7 @@ void FixEmitSurf::perform_task()
       // for CONSTANT or VARIABLE: task narget is fraction of its surf's area
       //   scale fraction by np or npcurrent (variable evaluation)
       // ninsert = rounded-down (ntarget + random number)
-      
+
       if (npmode == FLOW) ntarget = tasks[i].ntarget;
       else if (npmode == CONSTANT) ntarget = np * tasks[i].ntarget;
       else if (npmode == VARIABLE) ntarget = npcurrent * tasks[i].ntarget;
diff --git a/src/fix_emit_surf.h b/src/fix_emit_surf.h
index ebcaff322..b31a00a61 100644
--- a/src/fix_emit_surf.h
+++ b/src/fix_emit_surf.h
@@ -32,7 +32,6 @@ class FixEmitSurf : public FixEmit {
   FixEmitSurf(class SPARTA *, int, char **);
   ~FixEmitSurf();
   void init();
-  void setup();
 
   void grid_changed();
 
@@ -45,7 +44,7 @@ class FixEmitSurf : public FixEmit {
   int npmode,np;    // npmode = FLOW,CONSTANT,VARIABLE
   int npvar;
   char *npstr;
-  
+
   // copies of data from other classes
 
   int dimension,nspecies;
diff --git a/src/fix_grid_check.cpp b/src/fix_grid_check.cpp
index cbe263c87..bc4fe16cb 100644
--- a/src/fix_grid_check.cpp
+++ b/src/fix_grid_check.cpp
@@ -199,7 +199,7 @@ void FixGridCheck::end_of_step()
     // for split cell, also verify particle is in correct sub cell
     // expensive, so only do this check if requested
     // if cell volume = zero, error has already been flagged
-    
+
     if (!outside_check) continue;
     if (cells[icell].nsurf == 0) continue;
     if (cinfo[icell].volume == 0.0) continue;
@@ -209,11 +209,11 @@ void FixGridCheck::end_of_step()
 
     // check that particle is outside surfs
     // if no xcell found, cannot check
-    
+
     pflag = grid->point_outside_surfs(icell,xcell);
     if (!pflag) continue;
     pflag = grid->outside_surfs(icell,x,xcell);
-    
+
     if (!pflag) {
       if (outflag == ERROR) {
         char str[128];
@@ -228,7 +228,7 @@ void FixGridCheck::end_of_step()
     }
 
     // check that particle is in correct split subcell
-    
+
     if (cells[icell].nsplit <= 0) {
       int subcell;
       splitcell = sinfo[cells[icell].isplit].icell;
diff --git a/src/grid.cpp b/src/grid.cpp
index 09d0b3f99..5c511cae7 100644
--- a/src/grid.cpp
+++ b/src/grid.cpp
@@ -151,7 +151,7 @@ Grid::Grid(SPARTA *sparta) : Pointers(sparta)
 
   cut2d = NULL;
   cut3d = NULL;
-  
+
   // allocate hash for cell IDs
 
   hash = new MyHash();
diff --git a/src/grid_surf.cpp b/src/grid_surf.cpp
index b94ed7c14..a2efb8deb 100644
--- a/src/grid_surf.cpp
+++ b/src/grid_surf.cpp
@@ -1831,20 +1831,20 @@ int Grid::point_outside_surfs_explicit(int icell, double *x)
 
   double minsize = MIN(hi[0]-lo[0],hi[1]-lo[1]);
   double displace = EPSSURF * minsize;
-  
+
   if (dim == 2) {
     int npoint;
     double cpath[4],a[2],b[2];
     double *norm;
     Surf::Line *line;
-    
+
     Surf::Line *lines = surf->lines;
 
     for (int i = 0; i < nsurf; i++) {
       line = &lines[csurfs[i]];
       if (line->transparent) continue;
       norm = line->norm;
-      
+
       npoint = cut2d->clip_external(line->p1,line->p2,lo,hi,cpath);
       if (npoint < 2) continue;
 
@@ -1861,7 +1861,7 @@ int Grid::point_outside_surfs_explicit(int icell, double *x)
       x[2] = 0.0;
       return 1;
     }
-    
+
   } else {
     int npoint;
     double cpath[24],a[2],b[2];
@@ -1874,10 +1874,10 @@ int Grid::point_outside_surfs_explicit(int icell, double *x)
       tri = &tris[csurfs[i]];
       if (tri->transparent) continue;
       norm = tri->norm;
-      
+
       npoint = cut3d->clip_external(tri->p1,tri->p2,tri->p3,lo,hi,cpath);
       if (npoint < 3) continue;
-                    
+
       int face = cut3d->sameface_external(&cpath[0],&cpath[3],&cpath[6],lo,hi);
       if (face) {
         if (face == 1 and norm[0] < 0.0) continue;
@@ -1900,7 +1900,7 @@ int Grid::point_outside_surfs_explicit(int icell, double *x)
   // means entire cell is actually outside or inside, just touched by surfs
   // if outside, caller does not need to call outside_surfs()
   // if inside, caller can detect that its flow volume = zero
-  
+
   return 0;
 }