Commit ed8cc827 authored by sjplimp's avatar sjplimp Committed by GitHub
Browse files

Merge pull request #211 from akohlmey/add-respa-to-fix-flow-gauss

Add respa support to fix flow/gauss
parents 27dac024 144e6a80
Loading
Loading
Loading
Loading
+19 −1
Original line number Diff line number Diff line
@@ -63,7 +63,7 @@ applied by GD before computing a pressure drop or comparing it to
other methods, such as the pump method "(Zhu)"_#Zhu. The pressure
correction is discussed and described in "(Strong)"_#Strong.

NOTE: For a complete example including the considerations discussed
For a complete example including the considerations discussed
above, see the examples/USER/flow_gauss directory.

NOTE: Only the flux of the atoms in group-ID will be conserved. If the
@@ -93,6 +93,19 @@ work on the system must have {fix_modify energy yes} set as well. This
includes thermostat fixes and any constraints that hold the positions
of wall atoms fixed, such as "fix spring/self"_fix_spring_self.html.

If this fix is used in a simulation with the "rRESPA"_run_style.html
integrator, the applied acceleration must be computed and applied at the same
rRESPA level as the interactions between the flowing fluid and the obstacle.
The rRESPA level at which the acceleration is applied can be changed using
the "fix_modify"_fix_modify.html {respa} option discussed below. If the
flowing fluid and the obstacle interact through multiple interactions that are
computed at different rRESPA levels, then there must be a separate flow/gauss
fix for each level. For example, if the flowing fluid and obstacle interact
through pairwise and long-range Coulomb interactions, which are computed at
rRESPA levels 3 and 4, respectively, then there must be two separate
flow/gauss fixes, one that specifies {fix_modify respa 3} and one with
{fix_modify respa 4}.

:line

[Restart, fix_modify, output, run start/stop, minimize info:]
@@ -109,6 +122,11 @@ fix to subtract the work done from the
system's potential energy as part of "thermodynamic
output"_thermo_style.html.

The "fix_modify"_fix_modify.html {respa} option is supported by this
fix. This allows the user to set at which level of the "rRESPA"_run_style.html
integrator the fix computes and adds the external acceleration. Default is the
outermost level.

This fix computes a global scalar and a global 3-vector of forces,
which can be accessed by various "output
commands"_Section_howto.html#howto_15.  The scalar is the negative of the
+29 −2
Original line number Diff line number Diff line
@@ -27,6 +27,7 @@
#include "domain.h"
#include "error.h"
#include "citeme.h"
#include "respa.h"

using namespace LAMMPS_NS;
using namespace FixConst;
@@ -63,6 +64,8 @@ FixFlowGauss::FixFlowGauss(LAMMPS *lmp, int narg, char **arg) :
  extvector = 1;
  size_vector = 3;
  global_freq = 1;    //data available every timestep
  respa_level_support = 1;
  //default respa level=outermost level is set in init()

  dimension = domain->dimension;

@@ -109,9 +112,23 @@ int FixFlowGauss::setmask()
  int mask = 0;
  mask |= POST_FORCE;
  mask |= THERMO_ENERGY;
  mask |= POST_FORCE_RESPA;
  return mask;
}

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

void FixFlowGauss::init()
{
  //if respa level specified by fix_modify, then override default (outermost)
  //if specified level too high, set to max level
  if (strstr(update->integrate_style,"respa")) {
    ilevel_respa = ((Respa *) update->integrate)->nlevels-1;
    if (respa_level >= 0)
      ilevel_respa = MIN(respa_level,ilevel_respa);
  }
}

/* ----------------------------------------------------------------------
   setup is called after the initial evaluation of forces before a run, so we
   must remove the total force here too
@@ -127,6 +144,12 @@ void FixFlowGauss::setup(int vflag)
  if (mTot <= 0.0)
    error->all(FLERR,"Invalid group mass in fix flow/gauss");

  if (strstr(update->integrate_style,"respa")) {
    ((Respa *) update->integrate)->copy_flevel_f(ilevel_respa);
    post_force_respa(vflag,ilevel_respa,0);
    ((Respa *) update->integrate)->copy_f_flevel(ilevel_respa);
  }
  else
    post_force(vflag);
}

@@ -148,7 +171,6 @@ void FixFlowGauss::post_force(int vflag)
  int ii,jj;

  //find the total force on all atoms

  //initialize to zero
  double f_thisProc[3];
  for (ii=0; ii<3; ii++)
@@ -201,6 +223,11 @@ void FixFlowGauss::post_force(int vflag)

}

void FixFlowGauss::post_force_respa(int vflag, int ilevel, int iloop)
{
  if (ilevel == ilevel_respa) post_force(vflag);
}

/* ----------------------------------------------------------------------
   negative of work done by this fix
   This is only computed if requested, either with fix_modify energy yes, or with the energy keyword. Otherwise returns 0.
+5 −2
Original line number Diff line number Diff line
@@ -30,10 +30,12 @@ FixStyle(flow/gauss,FixFlowGauss)
    public:
      FixFlowGauss(class LAMMPS *, int, char **);
      int setmask();
      void init();
      void setup(int);
      void post_force(int);
      void post_force_respa(int, int, int);
      double compute_scalar();
      double compute_vector(int n);
      void post_force(int);
      void setup(int);

    protected:
      int dimension;
@@ -44,6 +46,7 @@ FixStyle(flow/gauss,FixFlowGauss)
      double pe_tot;      //total added energy
      double dt;          //timestep
      bool workflag;      //if calculate work done by fix
      int ilevel_respa;

    };