jianmu
/
lammps
forked from lijiext/lammps
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

new fix flow/gauss command 

git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@15508 f3b2605a-c512-4ea7-a41b-209d697bcdaa
This commit is contained in:
sjplimp 2016-08-27 22:18:45 +00:00
parent 646d5bb1b9
commit 90e6032f97
9 changed files with 1056 additions and 17 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

41
doc/html/Section_commands.html
View File

@ -778,12 +778,12 @@ USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT.</p>
<a class="reference internal" href="Section_start.html#start-3"><span class="std std-ref">LAMMPS is built with the appropriate package</span></a>.</p>
<table border="1" class="docutils">
<colgroup>
<col width="23%" />
<col width="14%" />
<col width="22%" />
<col width="15%" />
<col width="15%" />
<col width="15%" />
<col width="20%" />
<col width="13%" />
</colgroup>
<tbody valign="top">
<tr class="row-odd"><td><a class="reference internal" href="fix_adapt_fep.html"><span class="doc">adapt/fep</span></a></td>
@ -798,56 +798,63 @@ USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT.</p>
<td><a class="reference internal" href="fix_eos_cv.html"><span class="doc">eos/cv</span></a></td>
<td><a class="reference internal" href="fix_eos_table.html"><span class="doc">eos/table</span></a></td>
<td><a class="reference internal" href="fix_eos_table_rx.html"><span class="doc">eos/table/rx</span></a></td>
<td><a class="reference internal" href="fix_gle.html"><span class="doc">gle</span></a></td>
<td><a class="reference internal" href="fix_flow_gauss.html"><span class="doc">flow/gauss</span></a></td>
</tr>
<tr class="row-odd"><td><a class="reference internal" href="fix_imd.html"><span class="doc">imd</span></a></td>
<tr class="row-odd"><td><a class="reference internal" href="fix_gle.html"><span class="doc">gle</span></a></td>
<td><a class="reference internal" href="fix_imd.html"><span class="doc">imd</span></a></td>
<td><a class="reference internal" href="fix_ipi.html"><span class="doc">ipi</span></a></td>
<td><a class="reference internal" href="fix_langevin_drude.html"><span class="doc">langevin/drude</span></a></td>
<td><a class="reference internal" href="fix_langevin_eff.html"><span class="doc">langevin/eff</span></a></td>
<td><a class="reference internal" href="fix_lb_fluid.html"><span class="doc">lb/fluid</span></a></td>
<td><a class="reference internal" href="fix_lb_momentum.html"><span class="doc">lb/momentum</span></a></td>
</tr>
<tr class="row-even"><td><a class="reference internal" href="fix_lb_pc.html"><span class="doc">lb/pc</span></a></td>
<tr class="row-even"><td><a class="reference internal" href="fix_lb_momentum.html"><span class="doc">lb/momentum</span></a></td>
<td><a class="reference internal" href="fix_lb_pc.html"><span class="doc">lb/pc</span></a></td>
<td><a class="reference internal" href="fix_lb_rigid_pc_sphere.html"><span class="doc">lb/rigid/pc/sphere</span></a></td>
<td><a class="reference internal" href="fix_lb_viscous.html"><span class="doc">lb/viscous</span></a></td>
<td><a class="reference internal" href="fix_meso.html"><span class="doc">meso</span></a></td>
<td><a class="reference internal" href="fix_manifoldforce.html"><span class="doc">manifoldforce</span></a></td>
<td><a class="reference internal" href="fix_meso_stationary.html"><span class="doc">meso/stationary</span></a></td>
</tr>
<tr class="row-odd"><td><a class="reference internal" href="fix_nve_manifold_rattle.html"><span class="doc">nve/manifold/rattle</span></a></td>
<tr class="row-odd"><td><a class="reference internal" href="fix_meso_stationary.html"><span class="doc">meso/stationary</span></a></td>
<td><a class="reference internal" href="fix_nve_manifold_rattle.html"><span class="doc">nve/manifold/rattle</span></a></td>
<td><a class="reference internal" href="fix_nvt_manifold_rattle.html"><span class="doc">nvt/manifold/rattle</span></a></td>
<td><a class="reference internal" href="fix_nh_eff.html"><span class="doc">nph/eff</span></a></td>
<td><a class="reference internal" href="fix_nh_eff.html"><span class="doc">npt/eff</span></a></td>
<td><a class="reference internal" href="fix_nve_eff.html"><span class="doc">nve/eff</span></a></td>
<td><a class="reference internal" href="fix_nh_eff.html"><span class="doc">nvt/eff</span></a></td>
</tr>
<tr class="row-even"><td><a class="reference internal" href="fix_nvt_sllod_eff.html"><span class="doc">nvt/sllod/eff</span></a></td>
<tr class="row-even"><td><a class="reference internal" href="fix_nh_eff.html"><span class="doc">nvt/eff</span></a></td>
<td><a class="reference internal" href="fix_nvt_sllod_eff.html"><span class="doc">nvt/sllod/eff</span></a></td>
<td><a class="reference internal" href="fix_phonon.html"><span class="doc">phonon</span></a></td>
<td><a class="reference internal" href="fix_pimd.html"><span class="doc">pimd</span></a></td>
<td><a class="reference internal" href="fix_qbmsst.html"><span class="doc">qbmsst</span></a></td>
<td><a class="reference internal" href="fix_qeq_reax.html"><span class="doc">qeq/reax</span></a></td>
<td><a class="reference internal" href="fix_qmmm.html"><span class="doc">qmmm</span></a></td>
</tr>
<tr class="row-odd"><td><a class="reference internal" href="fix_qtb.html"><span class="doc">qtb</span></a></td>
<tr class="row-odd"><td><a class="reference internal" href="fix_qmmm.html"><span class="doc">qmmm</span></a></td>
<td><a class="reference internal" href="fix_qtb.html"><span class="doc">qtb</span></a></td>
<td><a class="reference internal" href="fix_reax_bonds.html"><span class="doc">reax/c/bonds</span></a></td>
<td><a class="reference internal" href="fix_reaxc_species.html"><span class="doc">reax/c/species</span></a></td>
<td><a class="reference internal" href="fix_rx.html"><span class="doc">rx</span></a></td>
<td><a class="reference internal" href="fix_saed_vtk.html"><span class="doc">saed/vtk</span></a></td>
<td><a class="reference internal" href="fix_shardlow.html"><span class="doc">shardlow</span></a></td>
</tr>
<tr class="row-even"><td><a class="reference internal" href="fix_smd.html"><span class="doc">smd</span></a></td>
<tr class="row-even"><td><a class="reference internal" href="fix_shardlow.html"><span class="doc">shardlow</span></a></td>
<td><a class="reference internal" href="fix_smd.html"><span class="doc">smd</span></a></td>
<td><a class="reference internal" href="fix_smd_adjust_dt.html"><span class="doc">smd/adjust/dt</span></a></td>
<td><a class="reference internal" href="fix_smd_integrate_tlsph.html"><span class="doc">smd/integrate/tlsph</span></a></td>
<td><a class="reference internal" href="fix_smd_integrate_ulsph.html"><span class="doc">smd/integrate/ulsph</span></a></td>
<td><a class="reference internal" href="fix_smd_move_triangulated_surface.html"><span class="doc">smd/move/triangulated/surface</span></a></td>
<td><a class="reference internal" href="fix_smd_setvel.html"><span class="doc">smd/setvel</span></a></td>
</tr>
<tr class="row-odd"><td><a class="reference internal" href="fix_smd_tlsph_reference_configuration.html"><span class="doc">smd/tlsph/reference/configuration</span></a></td>
<tr class="row-odd"><td><a class="reference internal" href="fix_smd_setvel.html"><span class="doc">smd/setvel</span></a></td>
<td><a class="reference internal" href="fix_smd_tlsph_reference_configuration.html"><span class="doc">smd/tlsph/reference/configuration</span></a></td>
<td><a class="reference internal" href="fix_smd_wall_surface.html"><span class="doc">smd/wall/surface</span></a></td>
<td><a class="reference internal" href="fix_temp_rescale_eff.html"><span class="doc">temp/rescale/eff</span></a></td>
<td><a class="reference internal" href="fix_ti_rs.html"><span class="doc">ti/rs</span></a></td>
<td><a class="reference internal" href="fix_ti_spring.html"><span class="doc">ti/spring</span></a></td>
<td><a class="reference internal" href="fix_ttm.html"><span class="doc">ttm/mod</span></a></td>
</tr>
<tr class="row-even"><td><a class="reference internal" href="fix_ttm.html"><span class="doc">ttm/mod</span></a></td>
<td>&nbsp;</td>
<td>&nbsp;</td>
<td>&nbsp;</td>
<td>&nbsp;</td>
<td>&nbsp;</td>
</tr>
</tbody>
</table>

321
doc/html/fix_flow_gauss.html Normal file
View File

@ -0,0 +1,321 @@
<!DOCTYPE html>
<!--[if IE 8]><html class="no-js lt-ie9" lang="en" > <![endif]-->
<!--[if gt IE 8]><!--> <html class="no-js" lang="en" > <!--<![endif]-->
<head>
<meta charset="utf-8">
<meta name="viewport" content="width=device-width, initial-scale=1.0">
<title>fix flow/gauss command &mdash; LAMMPS documentation</title>
<link rel="stylesheet" href="_static/css/theme.css" type="text/css" />
<link rel="stylesheet" href="_static/sphinxcontrib-images/LightBox2/lightbox2/css/lightbox.css" type="text/css" />
<link rel="top" title="LAMMPS documentation" href="index.html"/>
<script src="_static/js/modernizr.min.js"></script>
</head>
<body class="wy-body-for-nav" role="document">
<div class="wy-grid-for-nav">
<nav data-toggle="wy-nav-shift" class="wy-nav-side">
<div class="wy-side-nav-search">
<a href="Manual.html" class="icon icon-home"> LAMMPS
</a>
<div role="search">
<form id="rtd-search-form" class="wy-form" action="search.html" method="get">
<input type="text" name="q" placeholder="Search docs" />
<input type="hidden" name="check_keywords" value="yes" />
<input type="hidden" name="area" value="default" />
</form>
</div>
</div>
<div class="wy-menu wy-menu-vertical" data-spy="affix" role="navigation" aria-label="main navigation">
<ul>
<li class="toctree-l1"><a class="reference internal" href="Section_intro.html">1. Introduction</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_start.html">2. Getting Started</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_commands.html">3. Commands</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_packages.html">4. Packages</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_accelerate.html">5. Accelerating LAMMPS performance</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_howto.html">6. How-to discussions</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_example.html">7. Example problems</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_perf.html">8. Performance &amp; scalability</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_tools.html">9. Additional tools</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_modify.html">10. Modifying &amp; extending LAMMPS</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_python.html">11. Python interface to LAMMPS</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_errors.html">12. Errors</a></li>
<li class="toctree-l1"><a class="reference internal" href="Section_history.html">13. Future and history</a></li>
</ul>
</div>
&nbsp;
</nav>
<section data-toggle="wy-nav-shift" class="wy-nav-content-wrap">
<nav class="wy-nav-top" role="navigation" aria-label="top navigation">
<i data-toggle="wy-nav-top" class="fa fa-bars"></i>
<a href="Manual.html">LAMMPS</a>
</nav>
<div class="wy-nav-content">
<div class="rst-content">
<div role="navigation" aria-label="breadcrumbs navigation">
<ul class="wy-breadcrumbs">
<li><a href="Manual.html">Docs</a> &raquo;</li>
<li>fix flow/gauss command</li>
<li class="wy-breadcrumbs-aside">
<a href="http://lammps.sandia.gov">Website</a>
<a href="Section_commands.html#comm">Commands</a>
</li>
</ul>
<hr/>
</div>
<div role="main" class="document" itemscope="itemscope" itemtype="http://schema.org/Article">
<div itemprop="articleBody">
<div class="section" id="fix-flow-gauss-command">
<span id="index-0"></span><h1>fix flow/gauss command</h1>
<div class="section" id="syntax">
<h2>Syntax</h2>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">fix</span> <span class="n">ID</span> <span class="n">group</span><span class="o">-</span><span class="n">ID</span> <span class="n">flow</span><span class="o">/</span><span class="n">gauss</span> <span class="n">xflag</span> <span class="n">yflag</span> <span class="n">zflag</span> <span class="n">keyword</span>
</pre></div>
</div>
<ul class="simple">
<li>ID, group-ID are documented in <a class="reference internal" href="fix.html"><span class="doc">fix</span></a> command</li>
<li>flow/gauss = style name of this fix command</li>
<li>xflag,yflag,zflag = 0 or 1</li>
</ul>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="mi">0</span> <span class="o">=</span> <span class="n">do</span> <span class="ow">not</span> <span class="n">conserve</span> <span class="n">current</span> <span class="ow">in</span> <span class="n">this</span> <span class="n">dimension</span>
<span class="mi">1</span> <span class="o">=</span> <span class="n">conserve</span> <span class="n">current</span> <span class="ow">in</span> <span class="n">this</span> <span class="n">dimension</span>
</pre></div>
</div>
<ul class="simple">
<li>zero or more keyword/value pairs may be appended</li>
<li>keyword = <em>energy</em></li>
</ul>
<pre class="literal-block">
<em>energy</em> value = no or yes
no = do not compute work done by this fix
yes = compute work done by this fix
</pre>
</div>
<div class="section" id="examples">
<h2>Examples</h2>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">fix</span> <span class="n">GD</span> <span class="n">fluid</span> <span class="n">flow</span><span class="o">/</span><span class="n">gauss</span> <span class="mi">1</span> <span class="mi">0</span> <span class="mi">0</span>
<span class="n">fix</span> <span class="n">GD</span> <span class="n">fluid</span> <span class="n">flow</span><span class="o">/</span><span class="n">gauss</span> <span class="mi">1</span> <span class="mi">1</span> <span class="mi">1</span> <span class="n">energy</span> <span class="n">yes</span>
</pre></div>
</div>
</div>
<div class="section" id="description">
<h2>Description</h2>
<p>This fix implements the Gaussian dynamics (GD) method to simulate a system at
constant mass flux <a class="reference internal" href="#strong"><span class="std std-ref">(Strong)</span></a>. GD is a nonequilibrium molecular
dynamics simulation method that can be used to study fluid flows through
pores, pipes, and channels. In its original implementation GD was used to
compute the pressure required to achieve a fixed mass flux through an opening.
The flux can be conserved in any combination of the directions, x, y, or z,
using xflag,yflag,zflag. This fix does not initialize a net flux through
a system, it only conserves the center-of-mass momentum that is present
when the fix is declared in the input script. Use the <a class="reference internal" href="velocity.html"><span class="doc">velocity</span></a>
command to generate an initial center-of-mass momentum.</p>
<p>GD applies an external fluctuating gravitational field that acts as a driving
force to keep the system away from equilibrium. To maintain steady state, a
profile-unbiased thermostat must be implemented to dissipate the heat that is
added by the driving force. <a class="reference internal" href="compute_temp_profile.html"><span class="doc">Compute temp/profile</span></a>
can be used to implement a profile-unbiased thermostat.</p>
<p>A common use of this fix is to compute a pressure drop across a pipe, pore, or
membrane. The pressure profile can be computed in LAMMPS with <a class="reference internal" href="compute_stress_atom.html"><span class="doc">compute stress/atom</span></a> and <a class="reference internal" href="fix_ave_chunk.html"><span class="doc">fix ave/chunk</span></a>,
or with the hardy method in <a class="reference internal" href="fix_atc.html"><span class="doc">fix atc</span></a>. Note that the simple
<a class="reference internal" href="compute_stress_atom.html"><span class="doc">compute stress/atom</span></a> method is only accurate away
from inhomogeneities in the fluid, such as fixed wall atoms. Further, the
computed pressure profile must be corrected for the acceleration applied by
GD before computing a pressure drop or comparing it to other methods, such as
the pump method <a class="reference internal" href="#zhu"><span class="std std-ref">(Zhu)</span></a>. The pressure correction is discussed and
described in <a class="reference internal" href="#strong"><span class="std std-ref">(Strong)</span></a>.</p>
<div class="admonition note">
<p class="first admonition-title">Note</p>
<p class="last">For a complete example including the considerations discussed above, see
the examples/USER/flow_gauss directory.</p>
</div>
<div class="admonition note">
<p class="first admonition-title">Note</p>
<p class="last">Only the flux of the atoms in group-ID will be conserved. If the
velocities of the group-ID atoms are coupled to the velocities of other atoms
in the simulation, the flux will not be conserved. For example, in a
simulation with fluid atoms and harmonically constrained wall atoms, if a
single thermostat is applied to group <em>all</em>, the fluid atom velocities will be
coupled to the wall atom velocities, and the flux will not be conserved. This
issue can be avoided by thermostatting the fluid and wall groups separately.</p>
</div>
<dl class="docutils">
<dt>Adding an acceleration to atoms does work on the system. This added energy</dt>
<dd>can be optionally subtracted from the potential energy for the thermodynamic</dd>
</dl>
<p>output (see below) to check that the timestep is small enough to conserve
energy. Since the applied acceleration is fluctuating in time, the work cannot
be computed from a potential. As a result, computing the work is slightly more
computationally expensive than usual, so it is not performed by default. To
invoke the work calculation, use the <em>energy</em> keyword. The
<a class="reference internal" href="fix_modify.html"><span class="doc">fix_modify</span></a> <em>energy</em> option also invokes the work
calculation, and overrides an <em>energy no</em> setting here. If neither <em>energy yes</em>
or <em>fix_modify energy yes</em> are set, the global scalar computed by the fix
will return zero.</p>
<div class="admonition note">
<p class="first admonition-title">Note</p>
<p class="last">In order to check energy conservation, any other fixes that do work on
the system must have <em>fix_modify energy yes</em> set as well. This includes
thermostat fixes and any constraints that hold the positions of wall atoms
fixed, such as <a class="reference internal" href="fix_spring_self.html"><span class="doc">fix spring/self</span></a>.</p>
</div>
</div>
<hr class="docutils" />
<div class="section" id="restart-fix-modify-output-run-start-stop-minimize-info">
<h2>Restart, fix_modify, output, run start/stop, minimize info</h2>
<p>This fix is part of the USER-MISC package. It is only enabled if
LAMMPS was built with that package. See the <a class="reference internal" href="Section_start.html#start-3"><span class="std std-ref">Making LAMMPS</span></a> section for more info.</p>
<p>No information about this fix is written to <a class="reference internal" href="restart.html"><span class="doc">binary restart files</span></a>.</p>
<p>The <a class="reference internal" href="fix_modify.html"><span class="doc">fix_modify</span></a> <em>energy</em> option is supported by this
fix to subtract the work done from the
system&#8217;s potential energy as part of <a class="reference internal" href="thermo_style.html"><span class="doc">thermodynamic output</span></a>.</p>
<p>This fix computes a global scalar and a global 3-vector of forces,
which can be accessed by various <a class="reference internal" href="Section_howto.html#howto-15"><span class="std std-ref">output commands</span></a>. The scalar is the negative of the
work done on the system, see above discussion. The vector is the total force
that this fix applied to the group of atoms on the current timestep.
The scalar and vector values calculated by this fix are &#8220;extensive&#8221;.</p>
<p>No parameter of this fix can be used with the <em>start/stop</em> keywords of
the <a class="reference internal" href="run.html"><span class="doc">run</span></a> command.</p>
</div>
<div class="section" id="restrictions">
<h2>Restrictions</h2>
<blockquote>
<div>none</div></blockquote>
</div>
<div class="section" id="related-commands">
<h2>Related commands</h2>
<p><a class="reference internal" href="fix_addforce.html"><span class="doc">fix addforce</span></a>, <a class="reference internal" href="compute_temp_profile.html"><span class="doc">compute temp/profile</span></a>, <a class="reference internal" href="velocity.html"><span class="doc">velocity</span></a></p>
</div>
<div class="section" id="default">
<h2>Default</h2>
<p>The option default for the <em>energy</em> keyword is energy = no.</p>
<hr class="docutils" />
<p id="strong"><strong>(Strong)</strong> Strong and Eaves, J. Phys. Chem. Lett. 7, 1907 (2016).</p>
<p id="evans"><strong>(Evans)</strong> Evans and Morriss, Phys. Rev. Lett. 56, 2172 (1986).</p>
<p id="zhu"><strong>(Zhu)</strong> Zhu, Tajkhorshid, and Schulten, Biophys. J. 83, 154 (2002).</p>
<hr class="docutils" />
</div>
</div>
</div>
</div>
<footer>
<hr/>
<div role="contentinfo">
<p>
&copy; Copyright 2013 Sandia Corporation.
</p>
</div>
Built with <a href="http://sphinx-doc.org/">Sphinx</a> using a <a href="https://github.com/snide/sphinx_rtd_theme">theme</a> provided by <a href="https://readthedocs.org">Read the Docs</a>.
</footer>
</div>
</div>
</section>
</div>
<script type="text/javascript">
var DOCUMENTATION_OPTIONS = {
URL_ROOT:'./',
VERSION:'',
COLLAPSE_INDEX:false,
FILE_SUFFIX:'.html',
HAS_SOURCE: true
};
</script>
<script type="text/javascript" src="_static/jquery.js"></script>
<script type="text/javascript" src="_static/underscore.js"></script>
<script type="text/javascript" src="_static/doctools.js"></script>
<script type="text/javascript" src="https://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML"></script>
<script type="text/javascript" src="_static/sphinxcontrib-images/LightBox2/lightbox2/js/jquery-1.11.0.min.js"></script>
<script type="text/javascript" src="_static/sphinxcontrib-images/LightBox2/lightbox2/js/lightbox.min.js"></script>
<script type="text/javascript" src="_static/sphinxcontrib-images/LightBox2/lightbox2-customize/jquery-noconflict.js"></script>
<script type="text/javascript" src="_static/js/theme.js"></script>
<script type="text/javascript">
jQuery(function () {
SphinxRtdTheme.StickyNav.enable();
});
</script>
</body>
</html>

1
doc/src/Section_commands.txt
View File

@ -622,6 +622,7 @@ package"_Section_start.html#start_3.
"eos/cv"_fix_eos_cv.html,
"eos/table"_fix_eos_table.html,
"eos/table/rx"_fix_eos_table_rx.html,
"flow/gauss"_fix_flow_gauss.html,
"gle"_fix_gle.html,
"imd"_fix_imd.html,
"ipi"_fix_ipi.html,

137
doc/src/fix_flow_gauss.txt Normal file
View File

@ -0,0 +1,137 @@
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
:link(lws,http://lammps.sandia.gov)
:link(ld,Manual.html)
:link(lc,Section_commands.html#comm)
:line
fix flow/gauss command :h3
[Syntax:]
fix ID group-ID flow/gauss xflag yflag zflag keyword :pre
ID, group-ID are documented in "fix"_fix.html command :ulb,l
flow/gauss = style name of this fix command :l
xflag,yflag,zflag = 0 or 1 :l
0 = do not conserve current in this dimension
1 = conserve current in this dimension :pre
zero or more keyword/value pairs may be appended :l
keyword = {energy} :l
{energy} value = no or yes
no = do not compute work done by this fix
yes = compute work done by this fix :pre
:ule
[Examples:]
fix GD fluid flow/gauss 1 0 0
fix GD fluid flow/gauss 1 1 1 energy yes :pre
[Description:]
This fix implements the Gaussian dynamics (GD) method to simulate a system at
constant mass flux "(Strong)"_#Strong. GD is a nonequilibrium molecular
dynamics simulation method that can be used to study fluid flows through
pores, pipes, and channels. In its original implementation GD was used to
compute the pressure required to achieve a fixed mass flux through an opening.
The flux can be conserved in any combination of the directions, x, y, or z,
using xflag,yflag,zflag. This fix does not initialize a net flux through
a system, it only conserves the center-of-mass momentum that is present
when the fix is declared in the input script. Use the "velocity"_velocity.html
command to generate an initial center-of-mass momentum.
GD applies an external fluctuating gravitational field that acts as a driving
force to keep the system away from equilibrium. To maintain steady state, a
profile-unbiased thermostat must be implemented to dissipate the heat that is
added by the driving force. "Compute temp/profile"_compute_temp_profile.html
can be used to implement a profile-unbiased thermostat.
A common use of this fix is to compute a pressure drop across a pipe, pore, or
membrane. The pressure profile can be computed in LAMMPS with "compute
stress/atom"_compute_stress_atom.html and "fix ave/chunk"_fix_ave_chunk.html,
or with the hardy method in "fix atc"_fix_atc.html. Note that the simple
"compute stress/atom"_compute_stress_atom.html method is only accurate away
from inhomogeneities in the fluid, such as fixed wall atoms. Further, the
computed pressure profile must be corrected for the acceleration 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 above, see
the examples/USER/flow_gauss directory.
NOTE: Only the flux of the atoms in group-ID will be conserved. If the
velocities of the group-ID atoms are coupled to the velocities of other atoms
in the simulation, the flux will not be conserved. For example, in a
simulation with fluid atoms and harmonically constrained wall atoms, if a
single thermostat is applied to group {all}, the fluid atom velocities will be
coupled to the wall atom velocities, and the flux will not be conserved. This
issue can be avoided by thermostatting the fluid and wall groups separately.
Adding an acceleration to atoms does work on the system. This added energy
can be optionally subtracted from the potential energy for the thermodynamic
output (see below) to check that the timestep is small enough to conserve
energy. Since the applied acceleration is fluctuating in time, the work cannot
be computed from a potential. As a result, computing the work is slightly more
computationally expensive than usual, so it is not performed by default. To
invoke the work calculation, use the {energy} keyword. The
"fix_modify"_fix_modify.html {energy} option also invokes the work
calculation, and overrides an {energy no} setting here. If neither {energy yes}
or {fix_modify energy yes} are set, the global scalar computed by the fix
will return zero.
NOTE: In order to check energy conservation, any other fixes that do 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.
:line
[Restart, fix_modify, output, run start/stop, minimize info:]
This fix is part of the USER-MISC package. It is only enabled if
LAMMPS was built with that package. See the "Making
LAMMPS"_Section_start.html#start_3 section for more info.
No information about this fix is written to "binary restart
files"_restart.html.
The "fix_modify"_fix_modify.html {energy} option is supported by this
fix to subtract the work done from the
system's potential energy as part of "thermodynamic
output"_thermo_style.html.
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
work done on the system, see above discussion. The vector is the total force
that this fix applied to the group of atoms on the current timestep.
The scalar and vector values calculated by this fix are "extensive".
No parameter of this fix can be used with the {start/stop} keywords of
the "run"_run.html command.
[Restrictions:] none
[Related commands:]
"fix addforce"_fix_addforce.html, "compute temp/profile"_compute_temp_profile.html, "velocity"_velocity.html
[Default:]
The option default for the {energy} keyword is energy = no.
:line
:link(Strong)
[(Strong)] Strong and Eaves, J. Phys. Chem. Lett. 7, 1907 (2016).
:link(Evans)
[(Evans)] Evans and Morriss, Phys. Rev. Lett. 56, 2172 (1986).
:link(Zhu)
[(Zhu)] Zhu, Tajkhorshid, and Schulten, Biophys. J. 83, 154 (2002).
:line

45
examples/USER/flow_gauss/README Normal file
View File

@ -0,0 +1,45 @@
The input script in.GD is an example simulation using Gaussian dynamics (GD).
The simulation is of a simple 2d Lennard-Jones fluid flowing through a pipe.
For details see online LAMMPS documentation and
Strong and Eaves, J. Phys. Chem. Lett. 7(10) 2016, p. 1907.
Note that the run times and box size are chosen to allow a fast example run.
They are not adequate for a real simulation.
The script has the following parts:
1) initialize variables
These can be modified to customize the simulation. Note that if the
pipe dimensions L or d are changed, the geometry should be checked
by visualizing the coordinates in all.init.lammpstrj.
2) create box
3) set up potential
4) create atoms
5) set up profile-unbiased thermostat (PUT)
see Evans and Morriss, Phys. Rev. Lett. 56(20) 1986, p. 2172
By default, this uses boxes which contain on average 8 molecules.
6) equilibrate without GD
7) initialize the center-of-mass velocity and run to achieve steady-state
The system is initialized with a uniform velocity profile, which
relaxes over the course of the simulation.
8) collect data
The data is output in several files:
GD.out contains the force that GD applies, and the flux in the x- and
y- directions. The output Jx should be equal to the value of
J set in section 1, which is 0.1 by default.
x_profiles contains the velocity, density, and pressure profiles in
the x-direction. The pressure profile is given by
(-1/2V)*(c_spa[1] + c_spa[2]), where V is the volume of a
slice. The pressure profile is computed with IK1, see
Todd, Evans, and Davis, Phys. Rev. E 52(2) 1995, p. 1627.
Note that to compare with the pump method, or to
compute a pressure drop, you must correct this pressure
profile as described in Strong 2016 above.
Vy_profile is the velocity profile inside the pipe along the
y-direction, u_x(y).

258
examples/USER/flow_gauss/in.GD Executable file
View File

@ -0,0 +1,258 @@
#LAMMPS input script
#in.GD
#see README for details
###############################################################################
#initialize variables
clear
#frequency for outputting info (timesteps)
variable dump_rate equal 50
variable thermo_rate equal 10
#equilibration time (timesteps)
variable equil equal 1000
#stabilization time (timesteps to reach steady-state)
variable stabil equal 1000
#data collection time (timesteps)
variable run equal 2000
#length of pipe
variable L equal 30
#width of pipe
variable d equal 20
#flux (mass/sigma*tau)
variable J equal 0.1
#simulation box dimensions
variable Lx equal 100
variable Ly equal 40
#bulk fluid density
variable dens equal 0.8
#lattice spacing for wall atoms
variable aWall equal 1.0 #1.7472
#timestep
variable ts equal 0.001
#temperature
variable T equal 2.0
#thermostat damping constant
variable tdamp equal ${ts}*100
units lj
dimension 2
atom_style atomic
###############################################################################
#create box
#create lattice with the spacing aWall
variable rhoWall equal ${aWall}^(-2)
lattice sq ${rhoWall}
#modify input dimensions to be multiples of aWall
variable L1 equal round($L/${aWall})*${aWall}
variable d1 equal round($d/${aWall})*${aWall}
variable Ly1 equal round(${Ly}/${aWall})*${aWall}
variable Lx1 equal round(${Lx}/${aWall})*${aWall}
#create simulation box
variable lx2 equal ${Lx1}/2
variable ly2 equal ${Ly1}/2
region simbox block -${lx2} ${lx2} -${ly2} ${ly2} 0 0.1 units box
create_box 2 simbox
#####################################################################
#set up potential
mass 1 1.0 #fluid atoms
mass 2 1.0 #wall atoms
pair_style lj/cut 2.5
pair_modify shift yes
pair_coeff 1 1 1.0 1.0 2.5
pair_coeff 1 2 1.0 1.0 1.12246
pair_coeff 2 2 0.0 0.0 0.0
timestep ${ts}
#####################################################################
#create atoms
#create wall atoms everywhere
create_atoms 2 box
#define region which is "walled off"
variable dhalf equal ${d1}/2
variable Lhalf equal ${L1}/2
region walltop block -${Lhalf} ${Lhalf} ${dhalf} EDGE -0.1 0.1 &
units box
region wallbot block -${Lhalf} ${Lhalf} EDGE -${dhalf} -0.1 0.1 &
units box
region outsidewall union 2 walltop wallbot side out
#remove wall atoms outside wall region
group outside region outsidewall
delete_atoms group outside
#remove wall atoms that aren't on edge of wall region
variable x1 equal ${Lhalf}-${aWall}
variable y1 equal ${dhalf}+${aWall}
region insideTop block -${x1} ${x1} ${y1} EDGE -0.1 0.1 units box
region insideBot block -${x1} ${x1} EDGE -${y1} -0.1 0.1 units box
region insideWall union 2 insideTop insideBot
group insideWall region insideWall
delete_atoms group insideWall
#define new lattice, to give correct fluid density
#y lattice const must be a multiple of aWall
variable atrue equal ${dens}^(-1/2)
variable ay equal round(${atrue}/${aWall})*${aWall}
#choose x lattice const to give correct density
variable ax equal (${ay}*${dens})^(-1)
#change Lx to be multiple of ax
variable Lx1 equal round(${Lx}/${ax})*${ax}
variable lx2 equal ${Lx1}/2
change_box all x final -${lx2} ${lx2} units box
#define new lattice
lattice custom ${dens} &
a1 ${ax} 0.0 0.0 a2 0.0 ${ay} 0.0 a3 0.0 0.0 1.0 &
basis 0.0 0.0 0.0
#fill in rest of box with bulk particles
variable delta equal 0.001
variable Ldelt equal ${Lhalf}+${delta}
variable dDelt equal ${dhalf}-${delta}
region left block EDGE -${Ldelt} EDGE EDGE -0.1 0.1 units box
region right block ${Ldelt} EDGE EDGE EDGE -0.1 0.1 units box
region pipe block -${Ldelt} ${Ldelt} -${dDelt} ${dDelt} -0.1 0.1 &
units box
region bulk union 3 left pipe right
create_atoms 1 region bulk
group bulk type 1
group wall type 2
#remove atoms that are too close to wall
delete_atoms overlap 0.9 bulk wall
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
neigh_modify exclude group wall wall
velocity bulk create $T 78915 dist gaussian rot yes mom yes loop geom
#####################################################################
#set up PUT
#see Evans and Morriss, Phys. Rev. Lett. 56(20) 1986, p. 2172
#average number of particles per box, Evans and Morriss used 2.0
variable NperBox equal 8.0
#calculate box sizes
variable boxSide equal sqrt(${NperBox}/${dens})
variable nX equal round(lx/${boxSide})
variable nY equal round(ly/${boxSide})
variable dX equal lx/${nX}
variable dY equal ly/${nY}
#temperature of fluid (excluding wall)
compute myT bulk temp
#profile-unbiased temperature of fluid
compute myTp bulk temp/profile 1 1 0 xy ${nX} ${nY}
#thermo setup
thermo ${thermo_rate}
thermo_style custom step c_myT c_myTp etotal press
#dump initial configuration
dump 55 all custom 1 all.init.lammpstrj id type x y z vx vy vz
dump 56 wall custom 1 wall.init.lammpstrj id type x y z
dump_modify 55 sort id
dump_modify 56 sort id
run 0
undump 55
undump 56
#####################################################################
#equilibrate without GD
fix nvt bulk nvt temp $T $T ${tdamp}
fix_modify nvt temp myTp
fix 2 bulk enforce2d
run ${equil}
#####################################################################
#initialize the COM velocity and run to achieve steady-state
#calculate velocity to add: V=J/rho_total
variable Vadd equal $J*lx*ly/count(bulk)
#first remove any COM velocity, then add back the streaming velocity
velocity bulk zero linear
velocity bulk set ${Vadd} 0.0 0.0 units box sum yes mom no
fix GD bulk flow/gauss 1 0 0 #energy yes
#fix_modify GD energy yes
run ${stabil}
#####################################################################
#collect data
#print the applied force and total flux to ensure conservation of Jx
variable Fapp equal f_GD[1]
compute vxBulk bulk reduce sum vx
compute vyBulk bulk reduce sum vy
variable invVol equal 1.0/(lx*ly)
variable jx equal c_vxBulk*${invVol}
variable jy equal c_vyBulk*${invVol}
variable curr_step equal step
fix print_vCOM all print ${dump_rate} &
"${curr_step} ${Fapp} ${jx} ${jy}" file GD.out screen no &
title "timestep Fapp Jx Jy"
#compute IK1 pressure profile
#see Todd, Evans, and Davis, Phys. Rev. E 52(2) 1995, p. 1627
#use profile-unbiased temperature to remove the streaming velocity
#from the kinetic part of the pressure
compute spa bulk stress/atom myTp
#for the pressure profile, use the same grid as the PUT
compute chunkX bulk chunk/atom bin/1d x lower ${dX} units box
#output pressure profile and other profiles
#the pressure profile is (-1/2V)*(c_spa[1] + c_spa[2]), where
#V is the volume of a slice
fix profiles bulk ave/chunk 1 1 ${dump_rate} chunkX &
vx density/mass c_spa[1] c_spa[2] &
file x_profiles ave running overwrite
#compute velocity profile across the pipe with a finer grid
variable dYnew equal ${dY}/10
compute chunkY bulk chunk/atom bin/1d y center ${dYnew} units box &
region pipe
fix velYprof bulk ave/chunk 1 1 ${dump_rate} chunkY &
vx file Vy_profile ave running overwrite
#full trajectory
dump 7 bulk custom ${dump_rate} bulk.lammpstrj &
id type x y z
dump_modify 7 sort id
run ${run}

1
src/USER-MISC/README
View File

@ -37,6 +37,7 @@ dihedral_style spherical, Andrew Jewett, jewett.aij@gmail.com, 15 Jul 16
dihedral_style table, Andrew Jewett, jewett.aij@gmail.com, 10 Jan 12
fix addtorque, Laurent Joly (U Lyon), ljoly.ulyon at gmail.com, 8 Aug 11
fix ave/correlate/long, Jorge Ramirez (UPM Madrid), jorge.ramirez at upm.es, 21 Oct 2015
fix flow/gauss, Joel D. Eaves (CU Boulder), Joel.Eaves@Colorado.edu, 23 Aug 2016
fix gle, Michele Ceriotti (EPFL Lausanne), michele.ceriotti at gmail.com, 24 Nov 2014
fix imd, Axel Kohlmeyer, akohlmey at gmail.com, 9 Nov 2009
fix ipi, Michele Ceriotti (EPFL Lausanne), michele.ceriotti at gmail.com, 24 Nov 2014

216
src/USER-MISC/fix_flow_gauss.cpp Normal file
View File

@ -0,0 +1,216 @@
/* ----------------------------------------------------------------------
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.
Contributing authors: Steven E. Strong and Joel D. Eaves
Joel.Eaves@Colorado.edu
------------------------------------------------------------------------- */
#include <stdlib.h>
#include <string.h>
#include "fix_flow_gauss.h"
#include "atom.h"
#include "force.h"
#include "group.h"
#include "comm.h"
#include "update.h"
#include "domain.h"
#include "error.h"
#include "citeme.h"
using namespace LAMMPS_NS;
using namespace FixConst;
static const char cite_flow_gauss[] =
"Gaussian dynamics package:\n\n"
"@Article{strong_atomistic_2016,\n"
"title = {Atomistic Hydrodynamics and the Dynamical Hydrophobic Effect in Porous Graphene},\n"
"volume = {7},\n"
"number = {10},\n"
"issn = {1948-7185},\n"
"url = {http://dx.doi.org/10.1021/acs.jpclett.6b00748},\n"
"doi = {10.1021/acs.jpclett.6b00748},\n"
"urldate = {2016-05-10},\n"
"journal = {J. Phys. Chem. Lett.},\n"
"author = {Strong, Steven E. and Eaves, Joel D.},\n"
"year = {2016},\n"
"pages = {1907--1912}\n"
"}\n\n";
FixFlowGauss::FixFlowGauss(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
{
if (lmp->citeme) lmp->citeme->add(cite_flow_gauss);
if (narg < 6) error->all(FLERR,"Not enough input arguments");
// a group which conserves momentum must also conserve particle number
dynamic_group_allow = 0;
scalar_flag = 1;
vector_flag = 1;
extscalar = 1;
extvector = 1;
size_vector = 3;
global_freq = 1; //data available every timestep
dimension = domain->dimension;
//get inputs
int tmpFlag;
for (int ii=0; ii<3; ii++)
{
tmpFlag=force->inumeric(FLERR,arg[3+ii]);
if (tmpFlag==1 || tmpFlag==0)
flow[ii]=tmpFlag;
else
error->all(FLERR,"Constraint flags must be 1 or 0");
}
//by default, do not compute work done
workflag=0;
// process optional keyword
int iarg = 6;
while (iarg < narg) {
if ( strcmp(arg[iarg],"energy") == 0 ) {
if ( iarg+2 > narg ) error->all(FLERR,"Illegal energy keyword");
if ( strcmp(arg[iarg+1],"yes") == 0 ) workflag = 1;
else if ( strcmp(arg[iarg+1],"no") == 1 ) error->all(FLERR,"Illegal energy keyword");
iarg += 2;
} else error->all(FLERR,"Illegal fix flow/gauss command");
}
//error checking
if (dimension == 2) {
if (flow[2])
error->all(FLERR,"Can't constrain z flow in 2d simulation");
}
dt=update->dt;
pe_tot=0.0;
}
/* ---------------------------------------------------------------------- */
int FixFlowGauss::setmask()
{
int mask = 0;
mask |= POST_FORCE;
mask |= THERMO_ENERGY;
return mask;
}
/* ----------------------------------------------------------------------
setup is called after the initial evaluation of forces before a run, so we
must remove the total force here too
------------------------------------------------------------------------- */
void FixFlowGauss::setup(int vflag)
{
//need to compute work done if set fix_modify energy yes
if (thermo_energy)
workflag=1;
//get total mass of group
mTot=group->mass(igroup);
if (mTot <= 0.0)
error->all(FLERR,"Invalid group mass in fix flow/gauss");
post_force(vflag);
}
/* ----------------------------------------------------------------------
this is where Gaussian dynamics constraint is applied
------------------------------------------------------------------------- */
void FixFlowGauss::post_force(int vflag)
{
double **f = atom->f;
double **x = atom->x;
double **v = atom->v;
int *mask = atom->mask;
int *type = atom->type;
double *mass = atom->mass;
double *rmass = atom->rmass;
int nlocal = atom->nlocal;
int ii,jj;
//find the total force on all atoms
//initialize to zero
double f_thisProc[3];
for (ii=0; ii<3; ii++)
f_thisProc[ii]=0.0;
//add all forces on each processor
for(ii=0; ii<nlocal; ii++)
if (mask[ii] & groupbit)
for (jj=0; jj<3; jj++)
if (flow[jj])
f_thisProc[jj] += f[ii][jj];
//add the processor sums together
MPI_Allreduce(f_thisProc, f_tot, 3, MPI_DOUBLE, MPI_SUM, world);
//compute applied acceleration
for (ii=0; ii<3; ii++)
a_app[ii] = -f_tot[ii] / mTot;
//apply added accelleration to each atom
double f_app[3];
double peAdded=0.0;
for( ii = 0; ii<nlocal; ii++)
if (mask[ii] & groupbit) {
if (rmass) {
f_app[0] = a_app[0]*rmass[ii];
f_app[1] = a_app[1]*rmass[ii];
f_app[2] = a_app[2]*rmass[ii];
} else {
f_app[0] = a_app[0]*mass[type[ii]];
f_app[1] = a_app[1]*mass[type[ii]];
f_app[2] = a_app[2]*mass[type[ii]];
}
f[ii][0] += f_app[0]; //f_app[jj] is 0 if flow[jj] is false
f[ii][1] += f_app[1];
f[ii][2] += f_app[2];
//calculate added energy, since more costly, only do this if requested
if (workflag)
peAdded += f_app[0]*v[ii][0] + f_app[1]*v[ii][1] + f_app[2]*v[ii][2];
}
//finish calculation of work done, sum over all procs
if (workflag) {
double pe_tmp=0.0;
MPI_Allreduce(&peAdded,&pe_tmp,1,MPI_DOUBLE,MPI_SUM,world);
pe_tot += pe_tmp;
}
}
/* ----------------------------------------------------------------------
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.
------------------------------------------------------------------------- */
double FixFlowGauss::compute_scalar()
{
return -pe_tot*dt;
}
/* ----------------------------------------------------------------------
return components of applied force
------------------------------------------------------------------------- */
double FixFlowGauss::compute_vector(int n)
{
return -f_tot[n];
}

53
src/USER-MISC/fix_flow_gauss.h Normal file
View File

@ -0,0 +1,53 @@
/* -*- c++ -*- ----------------------------------------------------------
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.
Contributing authors: Steven E. Strong and Joel D. Eaves
Joel.Eaves@Colorado.edu
------------------------------------------------------------------------- */
#ifdef FIX_CLASS
FixStyle(flow/gauss,FixFlowGauss)
#else
#ifndef LMP_FIX_FLOWGAUSS_H
#define LMP_FIX_FLOWGAUSS_H
#include "fix.h"
namespace LAMMPS_NS {
class FixFlowGauss : public Fix {
public:
FixFlowGauss(class LAMMPS *, int, char **);
int setmask();
double compute_scalar();
double compute_vector(int n);
void post_force(int);
void setup(int);
protected:
int dimension;
bool flow[3]; //flag if each direction is conserved
double a_app[3]; //applied acceleration
double mTot; //total mass of constrained group
double f_tot[3]; //total applied force
double pe_tot; //total added energy
double dt; //timestep
bool workflag; //if calculate work done by fix
};
}
#endif
#endif