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

Merge pull request #27 from akohlmey/small-bugfixes 

Collected small changes and bugfixes
This commit is contained in:
sjplimp 2016-09-13 17:11:23 -06:00 committed by GitHub
parent dbc548dd88 80c5b01bfc
commit a69e059be3
19 changed files with 241 additions and 91 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

10
doc/html/boundary.html
View File

@ -173,7 +173,15 @@ you have set the <a class="reference internal" href="thermo_modify.html"><span c
allow for lost atoms.</p>
<p>For style <em>s</em>, the position of the face is set so as to encompass the
atoms in that dimension (shrink-wrapping), no matter how far they
move.</p>
move. Note that when the difference between the current box dimensions
and the shrink-wrap box dimensions is large, this can lead to lost
atoms at the beginning of a run when running in parallel. This is due
to the large change in the (global) box dimensions also causing
significant changes in the individual sub-domain sizes. If these
changes are farther than the communication cutoff, atoms will be lost.
This is best addressed by setting initial box dimensions to match the
shrink-wrapped dimensions more closely, by using <em>m</em> style boundaries
(see below).</p>
<p>For style <em>m</em>, shrink-wrapping occurs, but is bounded by the value
specified in the data or restart file or set by the
<a class="reference internal" href="create_box.html"><span class="doc">create_box</span></a> command. For example, if the upper z

24
doc/html/compute_rdf.html
View File

@ -141,13 +141,13 @@
</div>
<div class="section" id="examples">
<h2>Examples</h2>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">compute</span> <span class="mi">1</span> <span class="nb">all</span> <span class="n">rdf</span> <span class="mi">100</span>
<span class="n">compute</span> <span class="mi">1</span> <span class="nb">all</span> <span class="n">rdf</span> <span class="mi">100</span> <span class="mi">1</span> <span class="mi">1</span>
<span class="n">compute</span> <span class="mi">1</span> <span class="nb">all</span> <span class="n">rdf</span> <span class="mi">100</span> <span class="o">*</span> <span class="mi">3</span>
<span class="n">compute</span> <span class="mi">1</span> <span class="n">fluid</span> <span class="n">rdf</span> <span class="mi">500</span> <span class="mi">1</span> <span class="mi">1</span> <span class="mi">1</span> <span class="mi">2</span> <span class="mi">2</span> <span class="mi">1</span> <span class="mi">2</span> <span class="mi">2</span>
<span class="n">compute</span> <span class="mi">1</span> <span class="n">fluid</span> <span class="n">rdf</span> <span class="mi">500</span> <span class="mi">1</span><span class="o">*</span><span class="mi">3</span> <span class="mi">2</span> <span class="mi">5</span> <span class="o">*</span><span class="mi">10</span>
</pre></div>
</div>
<pre class="literal-block">
compute 1 all rdf 100
compute 1 all rdf 100 1 1
compute 1 all rdf 100 * 3
compute 1 fluid rdf 500 1 1 1 2 2 1 2 2
compute 1 fluid rdf 500 1*3 2 5 *10
</pre>
</div>
<div class="section" id="description">
<h2>Description</h2>
@ -184,7 +184,7 @@ listed, then a separate histogram is generated for each
<p>The <em>itypeN</em> and <em>jtypeN</em> settings can be specified in one of two
ways. An explicit numeric value can be used, as in the 4th example
above. Or a wild-card asterisk can be used to specify a range of atom
types. This takes the form &#8220;*&#8221; or &#8220;<em>n&#8221; or &#8220;n</em>&#8221; or &#8220;m*n&#8221;. If N = the
types. This takes the form &#8220;*&#8221; or &#8220;*n&#8221; or &#8220;n*&#8221; or &#8220;m*n&#8221;. If N = the
number of atom types, then an asterisk with no numeric values means
all types from 1 to N. A leading asterisk means all types from 1 to n
(inclusive). A trailing asterisk means all types from n to N
@ -226,10 +226,10 @@ atoms of type jtypeN.</p>
<p>The simplest way to output the results of the compute rdf calculation
to a file is to use the <a class="reference internal" href="fix_ave_time.html"><span class="doc">fix ave/time</span></a> command, for
example:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">compute</span> <span class="n">myRDF</span> <span class="nb">all</span> <span class="n">rdf</span> <span class="mi">50</span>
<span class="n">fix</span> <span class="mi">1</span> <span class="nb">all</span> <span class="n">ave</span><span class="o">/</span><span class="n">time</span> <span class="mi">100</span> <span class="mi">1</span> <span class="mi">100</span> <span class="n">c_myRDF</span><span class="p">[</span><span class="o">*</span><span class="p">]</span> <span class="n">file</span> <span class="n">tmp</span><span class="o">.</span><span class="n">rdf</span> <span class="n">mode</span> <span class="n">vector</span>
</pre></div>
</div>
<pre class="literal-block">
compute myRDF all rdf 50
fix 1 all ave/time 100 1 100 c_myRDF[*] file tmp.rdf mode vector
</pre>
<p><strong>Output info:</strong></p>
<p>This compute calculates a global array with the number of rows =
<em>Nbins</em>, and the number of columns = 1 + 2*Npairs, where Npairs is the

52
doc/html/fix_adapt.html
View File

@ -165,12 +165,12 @@
</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="mi">1</span> <span class="nb">all</span> <span class="n">adapt</span> <span class="mi">1</span> <span class="n">pair</span> <span class="n">soft</span> <span class="n">a</span> <span class="mi">1</span> <span class="mi">1</span> <span class="n">v_prefactor</span>
<span class="n">fix</span> <span class="mi">1</span> <span class="nb">all</span> <span class="n">adapt</span> <span class="mi">1</span> <span class="n">pair</span> <span class="n">soft</span> <span class="n">a</span> <span class="mi">2</span><span class="o">*</span> <span class="mi">3</span> <span class="n">v_prefactor</span>
<span class="n">fix</span> <span class="mi">1</span> <span class="nb">all</span> <span class="n">adapt</span> <span class="mi">1</span> <span class="n">pair</span> <span class="n">lj</span><span class="o">/</span><span class="n">cut</span> <span class="n">epsilon</span> <span class="o">*</span> <span class="o">*</span> <span class="n">v_scale1</span> <span class="n">coul</span><span class="o">/</span><span class="n">cut</span> <span class="n">scale</span> <span class="mi">3</span> <span class="mi">3</span> <span class="n">v_scale2</span> <span class="n">scale</span> <span class="n">yes</span> <span class="n">reset</span> <span class="n">yes</span>
<span class="n">fix</span> <span class="mi">1</span> <span class="nb">all</span> <span class="n">adapt</span> <span class="mi">10</span> <span class="n">atom</span> <span class="n">diameter</span> <span class="n">v_size</span>
</pre></div>
</div>
<pre class="literal-block">
fix 1 all adapt 1 pair soft a 1 1 v_prefactor
fix 1 all adapt 1 pair soft a 2* 3 v_prefactor
fix 1 all adapt 1 pair lj/cut epsilon * * v_scale1 coul/cut scale 3 3 v_scale2 scale yes reset yes
fix 1 all adapt 10 atom diameter v_size
</pre>
</div>
<div class="section" id="description">
<h2>Description</h2>
@ -225,8 +225,8 @@ meaning of these parameters:</p>
<table border="1" class="docutils">
<colgroup>
<col width="49%" />
<col width="35%" />
<col width="17%" />
<col width="36%" />
<col width="15%" />
</colgroup>
<tbody valign="top">
<tr class="row-odd"><td><a class="reference internal" href="pair_born.html"><span class="doc">born</span></a></td>
@ -257,23 +257,27 @@ meaning of these parameters:</p>
<td>epsilon,sigma,delta</td>
<td>type pairs</td>
</tr>
<tr class="row-even"><td><a class="reference internal" href="pair_lubricate.html"><span class="doc">lubricate</span></a></td>
<tr class="row-even"><td><a class="reference internal" href="pair_dipole.html"><span class="doc">lj/sf/dipole/sf</span></a></td>
<td>epsilon,sigma,scale</td>
<td>type pairs</td>
</tr>
<tr class="row-odd"><td><a class="reference internal" href="pair_lubricate.html"><span class="doc">lubricate</span></a></td>
<td>mu</td>
<td>global</td>
</tr>
<tr class="row-odd"><td><a class="reference internal" href="pair_gauss.html"><span class="doc">gauss</span></a></td>
<tr class="row-even"><td><a class="reference internal" href="pair_gauss.html"><span class="doc">gauss</span></a></td>
<td>a</td>
<td>type pairs</td>
</tr>
<tr class="row-even"><td><a class="reference internal" href="pair_morse.html"><span class="doc">morse</span></a></td>
<tr class="row-odd"><td><a class="reference internal" href="pair_morse.html"><span class="doc">morse</span></a></td>
<td>d0,r0,alpha</td>
<td>type pairs</td>
</tr>
<tr class="row-odd"><td><a class="reference internal" href="pair_soft.html"><span class="doc">soft</span></a></td>
<tr class="row-even"><td><a class="reference internal" href="pair_soft.html"><span class="doc">soft</span></a></td>
<td>a</td>
<td>type pairs</td>
</tr>
<tr class="row-even"><td><a class="reference internal" href="pair_kim.html"><span class="doc">kim</span></a></td>
<tr class="row-odd"><td><a class="reference internal" href="pair_kim.html"><span class="doc">kim</span></a></td>
<td>PARAM_FREE_*&amp;#58i,j,...</td>
<td>global</td>
</tr>
@ -315,7 +319,7 @@ each, as in the 1st example above. I &lt;= J is required. LAMMPS sets
the coefficients for the symmetric J,I interaction to the same values.</p>
<p>A wild-card asterisk can be used in place of or in conjunction with
the I,J arguments to set the coefficients for multiple pairs of atom
types. This takes the form &#8220;*&#8221; or &#8220;<em>n&#8221; or &#8220;n</em>&#8221; or &#8220;m*n&#8221;. If N = the
types. This takes the form &#8220;*&#8221; or &#8220;*n&#8221; or &#8220;n*&#8221; or &#8220;m*n&#8221;. If N = the
number of atom types, then an asterisk with no numeric values means
all types from 1 to N. A leading asterisk means all types from 1 to n
(inclusive). A trailing asterisk means all types from n to N
@ -341,10 +345,10 @@ details.</p>
<p>For example, these commands would change the prefactor coefficient of
the <a class="reference internal" href="pair_soft.html"><span class="doc">pair_style soft</span></a> potential from 10.0 to 30.0 in a
linear fashion over the course of a simulation:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">variable</span> <span class="n">prefactor</span> <span class="n">equal</span> <span class="n">ramp</span><span class="p">(</span><span class="mi">10</span><span class="p">,</span><span class="mi">30</span><span class="p">)</span>
<span class="n">fix</span> <span class="mi">1</span> <span class="nb">all</span> <span class="n">adapt</span> <span class="mi">1</span> <span class="n">pair</span> <span class="n">soft</span> <span class="n">a</span> <span class="o">*</span> <span class="o">*</span> <span class="n">v_prefactor</span>
</pre></div>
</div>
<pre class="literal-block">
variable prefactor equal ramp(10,30)
fix 1 all adapt 1 pair soft a * * v_prefactor
</pre>
<hr class="docutils" />
<p>The <em>kspace</em> keyword used the specified variable as a scale factor on
the energy, forces, virial calculated by whatever K-Space solver is
@ -380,14 +384,12 @@ constant).</p>
<p>For example, these commands would shrink the diameter of all granular
particles in the &#8220;center&#8221; group from 1.0 to 0.1 in a linear fashion
over the course of a 1000-step simulation:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">variable</span> <span class="n">size</span> <span class="n">equal</span> <span class="n">ramp</span><span class="p">(</span><span class="mf">1.0</span><span class="p">,</span><span class="mf">0.1</span><span class="p">)</span>
<span class="n">fix</span> <span class="mi">1</span> <span class="n">center</span> <span class="n">adapt</span> <span class="mi">10</span> <span class="n">atom</span> <span class="n">diameter</span> <span class="n">v_size</span>
</pre></div>
</div>
</div>
<pre class="literal-block">
variable size equal ramp(1.0,0.1)
fix 1 center adapt 10 atom diameter v_size
</pre>
<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><strong>Restart, fix_modify, output, run start/stop, minimize info:</strong></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>. None of the <a class="reference internal" href="fix_modify.html"><span class="doc">fix_modify</span></a> options
are relevant to this fix. No global or per-atom quantities are stored
by this fix for access by various <a class="reference internal" href="Section_howto.html#howto-15"><span class="std std-ref">output commands</span></a>. No parameter of this fix can

58
doc/html/pair_dipole.html
View File

@ -149,12 +149,12 @@
<h1>pair_style lj/long/dipole/long command</h1>
<div class="section" id="syntax">
<h2>Syntax</h2>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">pair_style</span> <span class="n">lj</span><span class="o">/</span><span class="n">cut</span><span class="o">/</span><span class="n">dipole</span><span class="o">/</span><span class="n">cut</span> <span class="n">cutoff</span> <span class="p">(</span><span class="n">cutoff2</span><span class="p">)</span>
<span class="n">pair_style</span> <span class="n">lj</span><span class="o">/</span><span class="n">sf</span><span class="o">/</span><span class="n">dipole</span><span class="o">/</span><span class="n">sf</span> <span class="n">cutoff</span> <span class="p">(</span><span class="n">cutoff2</span><span class="p">)</span>
<span class="n">pair_style</span> <span class="n">lj</span><span class="o">/</span><span class="n">cut</span><span class="o">/</span><span class="n">dipole</span><span class="o">/</span><span class="n">long</span> <span class="n">cutoff</span> <span class="p">(</span><span class="n">cutoff2</span><span class="p">)</span>
<span class="n">pair_style</span> <span class="n">lj</span><span class="o">/</span><span class="n">long</span><span class="o">/</span><span class="n">dipole</span><span class="o">/</span><span class="n">long</span> <span class="n">flag_lj</span> <span class="n">flag_coul</span> <span class="n">cutoff</span> <span class="p">(</span><span class="n">cutoff2</span><span class="p">)</span>
</pre></div>
</div>
<pre class="literal-block">
pair_style lj/cut/dipole/cut cutoff (cutoff2)
pair_style lj/sf/dipole/sf cutoff (cutoff2)
pair_style lj/cut/dipole/long cutoff (cutoff2)
pair_style lj/long/dipole/long flag_lj flag_coul cutoff (cutoff2)
</pre>
<ul class="simple">
<li>cutoff = global cutoff LJ (and Coulombic if only 1 arg) (distance units)</li>
<li>cutoff2 = global cutoff for Coulombic and dipole (optional) (distance units)</li>
@ -175,26 +175,27 @@
</div>
<div class="section" id="examples">
<h2>Examples</h2>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">pair_style</span> <span class="n">lj</span><span class="o">/</span><span class="n">cut</span><span class="o">/</span><span class="n">dipole</span><span class="o">/</span><span class="n">cut</span> <span class="mf">10.0</span>
<span class="n">pair_coeff</span> <span class="o">*</span> <span class="o">*</span> <span class="mf">1.0</span> <span class="mf">1.0</span>
<span class="n">pair_coeff</span> <span class="mi">2</span> <span class="mi">3</span> <span class="mf">1.0</span> <span class="mf">1.0</span> <span class="mf">2.5</span> <span class="mf">4.0</span>
</pre></div>
</div>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">pair_style</span> <span class="n">lj</span><span class="o">/</span><span class="n">sf</span><span class="o">/</span><span class="n">dipole</span><span class="o">/</span><span class="n">sf</span> <span class="mf">9.0</span>
<span class="n">pair_coeff</span> <span class="o">*</span> <span class="o">*</span> <span class="mf">1.0</span> <span class="mf">1.0</span>
<span class="n">pair_coeff</span> <span class="mi">2</span> <span class="mi">3</span> <span class="mf">1.0</span> <span class="mf">1.0</span> <span class="mf">2.5</span> <span class="mf">4.0</span>
</pre></div>
</div>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">pair_style</span> <span class="n">lj</span><span class="o">/</span><span class="n">cut</span><span class="o">/</span><span class="n">dipole</span><span class="o">/</span><span class="n">long</span> <span class="mf">10.0</span>
<span class="n">pair_coeff</span> <span class="o">*</span> <span class="o">*</span> <span class="mf">1.0</span> <span class="mf">1.0</span>
<span class="n">pair_coeff</span> <span class="mi">2</span> <span class="mi">3</span> <span class="mf">1.0</span> <span class="mf">1.0</span> <span class="mf">2.5</span> <span class="mf">4.0</span>
</pre></div>
</div>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">pair_style</span> <span class="n">lj</span><span class="o">/</span><span class="n">long</span><span class="o">/</span><span class="n">dipole</span><span class="o">/</span><span class="n">long</span> <span class="n">long</span> <span class="n">long</span> <span class="mf">3.5</span> <span class="mf">10.0</span>
<span class="n">pair_coeff</span> <span class="o">*</span> <span class="o">*</span> <span class="mf">1.0</span> <span class="mf">1.0</span>
<span class="n">pair_coeff</span> <span class="mi">2</span> <span class="mi">3</span> <span class="mf">1.0</span> <span class="mf">1.0</span> <span class="mf">2.5</span> <span class="mf">4.0</span>
</pre></div>
</div>
<pre class="literal-block">
pair_style lj/cut/dipole/cut 10.0
pair_coeff * * 1.0 1.0
pair_coeff 2 3 1.0 1.0 2.5 4.0
</pre>
<pre class="literal-block">
pair_style lj/sf/dipole/sf 9.0
pair_coeff * * 1.0 1.0
pair_coeff 2 3 1.0 1.0 2.5 4.0 scale 0.5
pair_coeff 2 3 1.0 1.0 2.5 4.0
</pre>
<pre class="literal-block">
pair_style lj/cut/dipole/long 10.0
pair_coeff * * 1.0 1.0
pair_coeff 2 3 1.0 1.0 2.5 4.0
</pre>
<pre class="literal-block">
pair_style lj/long/dipole/long long long 3.5 10.0
pair_coeff * * 1.0 1.0
pair_coeff 2 3 1.0 1.0 2.5 4.0
</pre>
</div>
<div class="section" id="description">
<h2>Description</h2>
@ -242,7 +243,10 @@ C.3 of <a class="reference internal" href="pair_gayberne.html#allen"><span class
<p>If one cutoff is specified in the pair_style command, it is used for
both the LJ and Coulombic (q,p) terms. If two cutoffs are specified,
they are used as cutoffs for the LJ and Coulombic (q,p) terms
respectively.</p>
respectively. This pair style also supports an optional <em>scale</em> keyword
as part of a pair_coeff statement, where the interactions can be
scaled according to this factor. This scale factor is also made available
for use with fix adapt.</p>
<p>Style <em>lj/cut/dipole/long</em> computes long-range point-dipole
interactions as discussed in <a class="reference internal" href="#toukmaji"><span class="std std-ref">(Toukmaji)</span></a>. Dipole-dipole,
dipole-charge, and charge-charge interactions are all supported, along

10
doc/src/boundary.txt
View File

@ -54,7 +54,15 @@ allow for lost atoms.
For style {s}, the position of the face is set so as to encompass the
atoms in that dimension (shrink-wrapping), no matter how far they
move.
move. Note that when the difference between the current box dimensions
and the shrink-wrap box dimensions is large, this can lead to lost
atoms at the beginning of a run when running in parallel. This is due
to the large change in the (global) box dimensions also causing
significant changes in the individual sub-domain sizes. If these
changes are farther than the communication cutoff, atoms will be lost.
This is best addressed by setting initial box dimensions to match the
shrink-wrapped dimensions more closely, by using {m} style boundaries
(see below).
For style {m}, shrink-wrapping occurs, but is bounded by the value
specified in the data or restart file or set by the

1
doc/src/fix_adapt.txt
View File

@ -110,6 +110,7 @@ meaning of these parameters:
"coul/long"_pair_coul.html: scale: type pairs:
"lj/cut"_pair_lj.html: epsilon,sigma: type pairs:
"lj/expand"_pair_lj_expand.html: epsilon,sigma,delta: type pairs:
"lj/sf/dipole/sf"_pair_dipole.html: epsilon,sigma,scale: type pairs:
"lubricate"_pair_lubricate.html: mu: global:
"gauss"_pair_gauss.html: a: type pairs:
"morse"_pair_morse.html: d0,r0,alpha: type pairs:

6
doc/src/pair_dipole.txt
View File

@ -41,6 +41,7 @@ pair_coeff 2 3 1.0 1.0 2.5 4.0 :pre
pair_style lj/sf/dipole/sf 9.0
pair_coeff * * 1.0 1.0
pair_coeff 2 3 1.0 1.0 2.5 4.0 scale 0.5
pair_coeff 2 3 1.0 1.0 2.5 4.0 :pre
pair_style lj/cut/dipole/long 10.0
@ -103,7 +104,10 @@ C.3 of "(Allen)"_#Allen.
If one cutoff is specified in the pair_style command, it is used for
both the LJ and Coulombic (q,p) terms. If two cutoffs are specified,
they are used as cutoffs for the LJ and Coulombic (q,p) terms
respectively.
respectively. This pair style also supports an optional {scale} keyword
as part of a pair_coeff statement, where the interactions can be
scaled according to this factor. This scale factor is also made available
for use with fix adapt.
Style {lj/cut/dipole/long} computes long-range point-dipole
interactions as discussed in "(Toukmaji)"_#Toukmaji. Dipole-dipole,

1
examples/COUPLE/fortran2/.gitignore vendored Normal file
View File

@ -0,0 +1 @@
*.mod

1
lib/.gitignore vendored Normal file
View File

@ -0,0 +1 @@
Makefile.lammps

1
lib/meam/.gitignore vendored Normal file
View File

@ -0,0 +1 @@
*.mod

1
src/MEAM/pair_meam.cpp
View File

@ -57,6 +57,7 @@ PairMEAM::PairMEAM(LAMMPS *lmp) : Pair(lmp)
arho1 = arho2 = arho3 = arho3b = t_ave = tsq_ave = NULL;
maxneigh = 0;
allocated = 0;
scrfcn = dscrfcn = fcpair = NULL;
nelements = 0;

2
src/RIGID/fix_rigid.cpp
View File

@ -2006,7 +2006,7 @@ void FixRigid::setup_bodies_static()
MPI_Allreduce(sum[0],all[0],6*nbody,MPI_DOUBLE,MPI_SUM,world);
// error check that re-computed momemts of inertia match diagonalized ones
// error check that re-computed moments of inertia match diagonalized ones
// do not do test for bodies with params read from infile
double norm;

2
src/RIGID/fix_rigid_small.cpp
View File

@ -2192,7 +2192,7 @@ void FixRigidSmall::setup_bodies_static()
commflag = ITENSOR;
comm->reverse_comm_fix(this,6);
// error check that re-computed momemts of inertia match diagonalized ones
// error check that re-computed moments of inertia match diagonalized ones
// do not do test for bodies with params read from infile
double norm;

2
src/USER-COLVARS/fix_colvars.cpp
View File

@ -383,7 +383,7 @@ void FixColvars::init()
error->all(FLERR,"Cannot use fix colvars without atom IDs");
if (atom->map_style == 0)
error->all(FLERR,"Fix colvars requires an atom map");
error->all(FLERR,"Fix colvars requires an atom map, see atom_modify");
if ((me == 0) && (update->whichflag == 2))
error->warning(FLERR,"Using fix colvars with minimization");

2
src/USER-COLVARS/fix_colvars.h
View File

@ -128,7 +128,7 @@ E: Cannot use fix colvars without atom IDs
Atom IDs are not defined, but fix colvars needs them to identify an atom.
E: Fix colvars requires an atom map
E: Fix colvars requires an atom map, see atom_modify
Use the atom_modify command to create an atom map.

127
src/USER-MISC/pair_lj_sf_dipole_sf.cpp
View File

@ -12,7 +12,8 @@
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Mario Orsi (U Southampton), orsimario@gmail.com
Contributing authors: Mario Orsi (QMUL), m.orsi@qmul.ac.uk
Samuel Genheden (University of Southampton)
------------------------------------------------------------------------- */
#include <math.h>
@ -34,7 +35,6 @@ using namespace LAMMPS_NS;
PairLJSFDipoleSF::PairLJSFDipoleSF(LAMMPS *lmp) : Pair(lmp)
{
single_enable = 0;
}
/* ---------------------------------------------------------------------- */
@ -55,6 +55,7 @@ PairLJSFDipoleSF::~PairLJSFDipoleSF()
memory->destroy(lj2);
memory->destroy(lj3);
memory->destroy(lj4);
memory->destroy(scale);
}
}
@ -86,6 +87,7 @@ void PairLJSFDipoleSF::compute(int eflag, int vflag)
double **torque = atom->torque;
int *type = atom->type;
int nlocal = atom->nlocal;
// The global scaling parameters aren't used anymore
double *special_coul = force->special_coul;
double *special_lj = force->special_lj;
int newton_pair = force->newton_pair;
@ -234,7 +236,7 @@ void PairLJSFDipoleSF::compute(int eflag, int vflag)
// total force
fq = factor_coul*qqrd2e;
fq = factor_coul*qqrd2e*scale[itype][jtype];
fx = fq*forcecoulx + delx*forcelj;
fy = fq*forcecouly + dely*forcelj;
fz = fq*forcecoulz + delz*forcelj;
@ -268,7 +270,7 @@ void PairLJSFDipoleSF::compute(int eflag, int vflag)
ecoul += -q[j] * r3inv * pqfac * pidotr;
if (mu[j][3] > 0.0 && qtmp != 0.0)
ecoul += qtmp * r3inv * qpfac * pjdotr;
ecoul *= factor_coul*qqrd2e;
ecoul *= factor_coul*qqrd2e*scale[itype][jtype];
} else ecoul = 0.0;
if (rsq < cut_ljsq[itype][jtype]) {
@ -315,6 +317,7 @@ void PairLJSFDipoleSF::allocate()
memory->create(lj2,n+1,n+1,"pair:lj2");
memory->create(lj3,n+1,n+1,"pair:lj3");
memory->create(lj4,n+1,n+1,"pair:lj4");
memory->create(scale,n+1,n+1,"pair:scale");
}
/* ----------------------------------------------------------------------
@ -352,7 +355,7 @@ void PairLJSFDipoleSF::settings(int narg, char **arg)
void PairLJSFDipoleSF::coeff(int narg, char **arg)
{
if (narg < 4 || narg > 6)
if (narg < 4 || narg > 7)
error->all(FLERR,"Incorrect args for pair coefficients");
if (!allocated) allocate();
@ -363,10 +366,13 @@ void PairLJSFDipoleSF::coeff(int narg, char **arg)
double epsilon_one = force->numeric(FLERR,arg[2]);
double sigma_one = force->numeric(FLERR,arg[3]);
double scale_one = 1.0;
if (narg >= 5) scale_one = force->numeric(FLERR,arg[4]);
double cut_lj_one = cut_lj_global;
double cut_coul_one = cut_coul_global;
if (narg >= 5) cut_coul_one = cut_lj_one = force->numeric(FLERR,arg[4]);
if (narg == 6) cut_coul_one = force->numeric(FLERR,arg[5]);
if (narg >= 6) cut_coul_one = cut_lj_one = force->numeric(FLERR,arg[5]);
if (narg == 7) cut_coul_one = force->numeric(FLERR,arg[6]);
int count = 0;
for (int i = ilo; i <= ihi; i++) {
@ -376,6 +382,7 @@ void PairLJSFDipoleSF::coeff(int narg, char **arg)
cut_lj[i][j] = cut_lj_one;
cut_coul[i][j] = cut_coul_one;
setflag[i][j] = 1;
scale[i][j] = scale_one;
count++;
}
}
@ -424,6 +431,7 @@ double PairLJSFDipoleSF::init_one(int i, int j)
lj2[j][i] = lj2[i][j];
lj3[j][i] = lj3[i][j];
lj4[j][i] = lj4[i][j];
scale[j][i] = scale[i][j];
return cut;
}
@ -445,6 +453,7 @@ void PairLJSFDipoleSF::write_restart(FILE *fp)
fwrite(&sigma[i][j],sizeof(double),1,fp);
fwrite(&cut_lj[i][j],sizeof(double),1,fp);
fwrite(&cut_coul[i][j],sizeof(double),1,fp);
fwrite(&scale[i][j],sizeof(double),1,fp);
}
}
}
@ -471,11 +480,13 @@ void PairLJSFDipoleSF::read_restart(FILE *fp)
fread(&sigma[i][j],sizeof(double),1,fp);
fread(&cut_lj[i][j],sizeof(double),1,fp);
fread(&cut_coul[i][j],sizeof(double),1,fp);
fread(&scale[i][j],sizeof(double),1,fp);
}
MPI_Bcast(&epsilon[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&sigma[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&cut_lj[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&cut_coul[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&scale[i][j],1,MPI_DOUBLE,0,world);
}
}
}
@ -506,3 +517,105 @@ void PairLJSFDipoleSF::read_restart_settings(FILE *fp)
MPI_Bcast(&cut_coul_global,1,MPI_DOUBLE,0,world);
MPI_Bcast(&mix_flag,1,MPI_INT,0,world);
}
// PairLJSFDipoleSF: calculation of force is missing (to be implemented)
double PairLJSFDipoleSF::single(int i, int j, int itype, int jtype, double rsq,
double factor_coul, double factor_lj,
double &fforce)
{
double r2inv,r6inv;
double pdotp,pidotr,pjdotr,pre1,delx,dely,delz;
double rinv, r3inv,r5inv, rcutlj2inv, rcutcoul2inv,rcutlj6inv;
double qtmp,xtmp,ytmp,ztmp,bfac,pqfac,qpfac, ecoul, evdwl;
double **x = atom->x;
double *q = atom->q;
double **mu = atom->mu;
qtmp = q[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
r2inv = 1.0/rsq;
rinv = sqrt(r2inv);
fforce = 0.0;
if (rsq < cut_coulsq[itype][jtype]) {
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
// if (qtmp != 0.0 && q[j] != 0.0) {
// pre1 = qtmp*q[j]*rinv*(r2inv-1.0/cut_coulsq[itype][jtype]);
// forcecoulx += pre1*delx;
// forcecouly += pre1*dely;
// forcecoulz += pre1*delz;
// }
if (mu[i][3] > 0.0 && mu[j][3] > 0.0) {
r3inv = r2inv*rinv;
r5inv = r3inv*r2inv;
rcutcoul2inv=1.0/cut_coulsq[itype][jtype];
pdotp = mu[i][0]*mu[j][0] + mu[i][1]*mu[j][1] + mu[i][2]*mu[j][2];
pidotr = mu[i][0]*delx + mu[i][1]*dely + mu[i][2]*delz;
pjdotr = mu[j][0]*delx + mu[j][1]*dely + mu[j][2]*delz;
bfac = 1.0 - 4.0*rsq*sqrt(rsq)*rcutcoul2inv*sqrt(rcutcoul2inv) +
3.0*rsq*rsq*rcutcoul2inv*rcutcoul2inv;
}
if (mu[i][3] > 0.0 && q[j] != 0.0) {
r3inv = r2inv*rinv;
r5inv = r3inv*r2inv;
pidotr = mu[i][0]*delx + mu[i][1]*dely + mu[i][2]*delz;
rcutcoul2inv=1.0/cut_coulsq[itype][jtype];
pqfac = 1.0 - 3.0*rsq*rcutcoul2inv +
2.0*rsq*sqrt(rsq)*rcutcoul2inv*sqrt(rcutcoul2inv);
}
if (mu[j][3] > 0.0 && qtmp != 0.0) {
r3inv = r2inv*rinv;
r5inv = r3inv*r2inv;
pjdotr = mu[j][0]*delx + mu[j][1]*dely + mu[j][2]*delz;
rcutcoul2inv=1.0/cut_coulsq[itype][jtype];
qpfac = 1.0 - 3.0*rsq*rcutcoul2inv +
2.0*rsq*sqrt(rsq)*rcutcoul2inv*sqrt(rcutcoul2inv);
}
}
if (rsq < cut_ljsq[itype][jtype]) {
r6inv = r2inv*r2inv*r2inv;
rcutlj2inv = 1.0 / cut_ljsq[itype][jtype];
rcutlj6inv = rcutlj2inv * rcutlj2inv * rcutlj2inv;
}
double eng = 0.0;
if (rsq < cut_coulsq[itype][jtype]) {
ecoul = (1.0-sqrt(rsq)/sqrt(cut_coulsq[itype][jtype]));
ecoul *= ecoul;
ecoul *= qtmp * q[j] * rinv;
if (mu[i][3] > 0.0 && mu[j][3] > 0.0)
ecoul += bfac * (r3inv*pdotp - 3.0*r5inv*pidotr*pjdotr);
if (mu[i][3] > 0.0 && q[j] != 0.0)
ecoul += -q[j] * r3inv * pqfac * pidotr;
if (mu[j][3] > 0.0 && qtmp != 0.0)
ecoul += qtmp * r3inv * qpfac * pjdotr;
ecoul *= factor_coul*force->qqrd2e*scale[itype][jtype];
eng += ecoul;
}
if (rsq < cut_ljsq[itype][jtype]) {
evdwl = r6inv*(lj3[itype][jtype]*r6inv-lj4[itype][jtype])+
rcutlj6inv*(6*lj3[itype][jtype]*rcutlj6inv-3*lj4[itype][jtype])*
rsq*rcutlj2inv+
rcutlj6inv*(-7*lj3[itype][jtype]*rcutlj6inv+4*lj4[itype][jtype]);
eng += evdwl*factor_lj;
}
return eng;
}
/* ---------------------------------------------------------------------- */
void *PairLJSFDipoleSF::extract(const char *str, int &dim)
{
dim = 2;
if (strcmp(str,"epsilon") == 0) return (void *) epsilon;
if (strcmp(str,"sigma") == 0) return (void *) sigma;
if (strcmp(str,"scale") == 0) return (void *) scale;
return NULL;
}

3
src/USER-MISC/pair_lj_sf_dipole_sf.h
View File

@ -37,6 +37,8 @@ class PairLJSFDipoleSF : public Pair {
void read_restart(FILE *);
void write_restart_settings(FILE *);
void read_restart_settings(FILE *);
virtual double single(int, int, int, int, double, double, double, double &);
void *extract(const char *, int &);
protected:
double cut_lj_global,cut_coul_global;
@ -44,6 +46,7 @@ class PairLJSFDipoleSF : public Pair {
double **cut_coul,**cut_coulsq;
double **epsilon,**sigma;
double **lj1,**lj2,**lj3,**lj4;
double **scale;
void allocate();
};

26
src/USER-PHONON/fix_phonon.cpp
View File

@ -190,7 +190,7 @@ FixPhonon::FixPhonon(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
sprintf(str,"Can not open output file %s",logfile);
error->one(FLERR,str);
}
for (int i = 0; i < 60; ++i) fprintf(flog,"#"); fprintf(flog,"\n");
fprintf(flog,"############################################################\n");
fprintf(flog,"# group name of the atoms under study : %s\n", group->names[igroup]);
fprintf(flog,"# total number of atoms in the group : %d\n", ngroup);
fprintf(flog,"# dimension of the system : %d D\n", sysdim);
@ -200,7 +200,7 @@ FixPhonon::FixPhonon(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
fprintf(flog,"# frequency of the measurement : %d\n", nevery);
fprintf(flog,"# output result after this many measurement: %d\n", nfreq);
fprintf(flog,"# number of processors used by this run : %d\n", nprocs);
for (int i = 0; i < 60; ++i) fprintf(flog,"#"); fprintf(flog,"\n");
fprintf(flog,"############################################################\n");
fprintf(flog,"# mapping information between lattice indices and atom id\n");
fprintf(flog,"# nx ny nz nucell\n");
fprintf(flog,"%d %d %d %d\n", nx, ny, nz, nucell);
@ -214,7 +214,7 @@ FixPhonon::FixPhonon(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
ix = (idx/(nucell*nz*ny))%nx;
fprintf(flog,"%d %d %d %d " TAGINT_FORMAT "\n", ix, iy, iz, iu, itag);
}
for (int i = 0; i < 60; ++i) fprintf(flog,"#"); fprintf(flog,"\n");
fprintf(flog,"############################################################\n");
fflush(flog);
}
surf2tag.clear();
@ -734,16 +734,16 @@ void FixPhonon::postprocess( )
fclose(fp_bin);
// write log file, here however, it is the dynamical matrix that is written
for (int i = 0; i < 60; ++i) fprintf(flog,"#"); fprintf(flog,"\n");
fprintf(flog, "# Current time step : " BIGINT_FORMAT "\n", update->ntimestep);
fprintf(flog, "# Total number of measurements : %d\n", neval);
fprintf(flog, "# Average temperature of the measurement : %lg\n", TempAve);
fprintf(flog, "# Boltzmann constant under current units : %lg\n", boltz);
fprintf(flog, "# basis vector A1 = [%lg %lg %lg]\n", basevec[0], basevec[1], basevec[2]);
fprintf(flog, "# basis vector A2 = [%lg %lg %lg]\n", basevec[3], basevec[4], basevec[5]);
fprintf(flog, "# basis vector A3 = [%lg %lg %lg]\n", basevec[6], basevec[7], basevec[8]);
for (int i = 0; i < 60; ++i) fprintf(flog,"#"); fprintf(flog,"\n");
fprintf(flog, "# qx\t qy \t qz \t\t Phi(q)\n");
fprintf(flog,"############################################################\n");
fprintf(flog,"# Current time step : " BIGINT_FORMAT "\n", update->ntimestep);
fprintf(flog,"# Total number of measurements : %d\n", neval);
fprintf(flog,"# Average temperature of the measurement : %lg\n", TempAve);
fprintf(flog,"# Boltzmann constant under current units : %lg\n", boltz);
fprintf(flog,"# basis vector A1 = [%lg %lg %lg]\n", basevec[0], basevec[1], basevec[2]);
fprintf(flog,"# basis vector A2 = [%lg %lg %lg]\n", basevec[3], basevec[4], basevec[5]);
fprintf(flog,"# basis vector A3 = [%lg %lg %lg]\n", basevec[6], basevec[7], basevec[8]);
fprintf(flog,"############################################################\n");
fprintf(flog,"# qx\t qy \t qz \t\t Phi(q)\n");
EnforceASR();

3
src/VORONOI/compute_voronoi_atom.cpp
View File

@ -127,6 +127,9 @@ ComputeVoronoi::ComputeVoronoi(LAMMPS *lmp, int narg, char **arg) :
if (occupation && ( surface!=VOROSURF_NONE || maxedge>0 ) )
error->all(FLERR,"Illegal compute voronoi/atom command (occupation and (surface or edges))");
if (occupation && (atom->map_style == 0))
error->all(FLERR,"Compute voronoi/atom occupation requires an atom map, see atom_modify");
nmax = rmax = 0;
edge = rfield = sendvector = NULL;
voro = NULL;