Merge branch 'develop' into collected-small-changes
This commit is contained in:
Axel Kohlmeyer
@ -1,7 +1,7 @@
|
||||
Use chunks to calculate system properties
|
||||
=========================================
|
||||
|
||||
In LAMMS, "chunks" are collections of atoms, as defined by the
|
||||
In LAMMPS, "chunks" are collections of atoms, as defined by the
|
||||
:doc:`compute chunk/atom <compute_chunk_atom>` command, which assigns
|
||||
each atom to a chunk ID (or to no chunk at all). The number of chunks
|
||||
and the assignment of chunk IDs to atoms can be static or change over
|
||||
@ -148,14 +148,14 @@ Example calculations with chunks
|
||||
Here are examples using chunk commands to calculate various
|
||||
properties:
|
||||
|
||||
(1) Average velocity in each of 1000 2d spatial bins:
|
||||
1. Average velocity in each of 1000 2d spatial bins:
|
||||
|
||||
.. code-block:: LAMMPS
|
||||
|
||||
compute cc1 all chunk/atom bin/2d x 0.0 0.1 y lower 0.01 units reduced
|
||||
fix 1 all ave/chunk 100 10 1000 cc1 vx vy file tmp.out
|
||||
|
||||
(2) Temperature in each spatial bin, after subtracting a flow
|
||||
2. Temperature in each spatial bin, after subtracting a flow
|
||||
velocity:
|
||||
|
||||
.. code-block:: LAMMPS
|
||||
@ -164,7 +164,7 @@ velocity:
|
||||
compute vbias all temp/profile 1 0 0 y 10
|
||||
fix 1 all ave/chunk 100 10 1000 cc1 temp bias vbias file tmp.out
|
||||
|
||||
(3) Center of mass of each molecule:
|
||||
3. Center of mass of each molecule:
|
||||
|
||||
.. code-block:: LAMMPS
|
||||
|
||||
@ -172,7 +172,7 @@ velocity:
|
||||
compute myChunk all com/chunk cc1
|
||||
fix 1 all ave/time 100 1 100 c_myChunk[*] file tmp.out mode vector
|
||||
|
||||
(4) Total force on each molecule and ave/max across all molecules:
|
||||
4. Total force on each molecule and ave/max across all molecules:
|
||||
|
||||
.. code-block:: LAMMPS
|
||||
|
||||
@ -183,7 +183,7 @@ velocity:
|
||||
thermo 1000
|
||||
thermo_style custom step temp v_xave v_xmax
|
||||
|
||||
(5) Histogram of cluster sizes:
|
||||
5. Histogram of cluster sizes:
|
||||
|
||||
.. code-block:: LAMMPS
|
||||
|
||||
@ -192,16 +192,16 @@ velocity:
|
||||
compute size all property/chunk cc1 count
|
||||
fix 1 all ave/histo 100 1 100 0 20 20 c_size mode vector ave running beyond ignore file tmp.histo
|
||||
|
||||
(6) An example for using a per-chunk value to apply per-atom forces to
|
||||
6. An example for using a per-chunk value to apply per-atom forces to
|
||||
compress individual polymer chains (molecules) in a mixture, is
|
||||
explained on the :doc:`compute chunk/spread/atom <compute_chunk_spread_atom>` command doc page.
|
||||
|
||||
(7) An example for using one set of per-chunk values for molecule
|
||||
7. An example for using one set of per-chunk values for molecule
|
||||
chunks, to create a second set of micelle-scale chunks (clustered
|
||||
molecules, due to hydrophobicity), is explained on the
|
||||
:doc:`compute reduce/chunk <compute_reduce_chunk>` command doc page.
|
||||
|
||||
(8) An example for using one set of per-chunk values (dipole moment
|
||||
8. An example for using one set of per-chunk values (dipole moment
|
||||
vectors) for molecule chunks, spreading the values to each atom in
|
||||
each chunk, then defining a second set of chunks as spatial bins, and
|
||||
using the :doc:`fix ave/chunk <fix_ave_chunk>` command to calculate an
|
||||
|
||||
@ -38,7 +38,7 @@ Syntax
|
||||
*electrode/thermo* args = potential eta *temp* values
|
||||
potential = electrode potential
|
||||
charge = electrode charge
|
||||
eta = reciprocal width of electrode charge smearing
|
||||
eta = reciprocal width of electrode charge smearing (can be NULL if eta keyword is used)
|
||||
*temp* values = T_v tau_v rng_v
|
||||
T_v = temperature of thermo-potentiostat
|
||||
tau_v = time constant of thermo-potentiostat
|
||||
@ -110,7 +110,7 @@ electrostatic configurations:
|
||||
:ref:`(Deissenbeck)<Deissenbeck>` between two electrodes
|
||||
|
||||
* (resulting in changing charges and potentials with appropriate
|
||||
average potential difference and thermal variance)
|
||||
average potential difference and thermal variance)
|
||||
|
||||
The first group-ID provided to each fix specifies the first electrode
|
||||
group, and more group(s) are added using the *couple* keyword for each
|
||||
@ -287,8 +287,18 @@ The *fix_modify tf* option enables the Thomas-Fermi metallicity model
|
||||
fix_modify ID tf type length voronoi
|
||||
|
||||
|
||||
If this option is used parameters must be set for all atom types of the
|
||||
electrode.
|
||||
If this option is used, these two parameters must be set for
|
||||
all atom types of the electrode:
|
||||
|
||||
* `tf` is the Thomas-Fermi length :math:`l_{TF}`
|
||||
* `voronoi` is the Voronoi volume per atom in units of length cubed
|
||||
|
||||
Different types may have different `tf` and `voronoi` values.
|
||||
The following self-energy term is then added for all electrode atoms:
|
||||
|
||||
.. math::
|
||||
|
||||
A_{ii} += \frac{1}{4 \pi \epsilon_0} \times \frac{4 \pi l_{TF}^2}{\mathrm{Voronoi volume}}
|
||||
|
||||
The *fix_modify timer* option turns on (off) additional timer outputs in the log
|
||||
file, for code developers to track optimization.
|
||||
@ -321,9 +331,11 @@ The global array has *N* rows and *2N+1* columns, where the fix manages
|
||||
array, the elements are:
|
||||
|
||||
* array[I][1] = total charge that group *I* would have had *if it were
|
||||
at 0 V applied potential* * array[I][2 to *N* + 1] = the *N* entries
|
||||
at 0 V applied potential*
|
||||
* array[I][2 to *N* + 1] = the *N* entries
|
||||
of the *I*-th row of the electrode capacitance matrix (definition
|
||||
follows) * array[I][*N* + 2 to *2N* + 1] = the *N* entries of the
|
||||
follows)
|
||||
* array[I][*N* + 2 to *2N* + 1] = the *N* entries of the
|
||||
*I*-th row of the electrode elastance matrix (the inverse of the
|
||||
electrode capacitance matrix)
|
||||
|
||||
|
||||
@ -12,7 +12,7 @@ Syntax
|
||||
|
||||
* file = name of data file to read in
|
||||
* zero or more keyword/arg pairs may be appended
|
||||
* keyword = *add* or *offset* or *shift* or *extra/atom/types* or *extra/bond/types* or *extra/angle/types* or *extra/dihedral/types* or *extra/improper/types* or *extra/bond/per/atom* or *extra/angle/per/atom* or *extra/dihedral/per/atom* or *extra/improper/per/atom* or *group* or *nocoeff* or *fix*
|
||||
* keyword = *add* or *offset* or *shift* or *extra/atom/types* or *extra/bond/types* or *extra/angle/types* or *extra/dihedral/types* or *extra/improper/types* or *extra/bond/per/atom* or *extra/angle/per/atom* or *extra/dihedral/per/atom* or *extra/improper/per/atom* or *extra/special/per/atom* or *group* or *nocoeff* or *fix*
|
||||
|
||||
.. parsed-literal::
|
||||
|
||||
|
||||
@ -32,7 +32,7 @@ Syntax
|
||||
|
||||
.. code-block:: LAMMPS
|
||||
|
||||
reset atoms mol group-ID keyword value ...
|
||||
reset_atoms mol group-ID keyword value ...
|
||||
|
||||
* group-ID = ID of group of atoms whose molecule IDs will be reset
|
||||
* zero or more keyword/value pairs can be appended
|
||||
@ -66,16 +66,16 @@ Description
|
||||
.. versionadded:: 22Dec2022
|
||||
|
||||
The *reset_atoms* command resets the values of a specified atom
|
||||
property. In contrast to the set command, it does this in a
|
||||
property. In contrast to the *set* command, it does this in a
|
||||
collective manner which resets the values for many atoms in a
|
||||
self-consistent way. This is often useful when the simulated system
|
||||
has undergone significant modifications like adding or removing atoms
|
||||
or molecules, joining data files, changing bonds, or large-scale
|
||||
self-consistent way. This command is often useful when the simulated
|
||||
system has undergone significant modifications like adding or removing
|
||||
atoms or molecules, joining data files, changing bonds, or large-scale
|
||||
diffusion.
|
||||
|
||||
The new values can be thought of as a *reset*, similar to values atoms
|
||||
would have if a new data file were being read or a new simulation
|
||||
performed. Note that the set command also resets atom properties to
|
||||
performed. Note that the *set* command also resets atom properties to
|
||||
new values, but it treats each atom independently.
|
||||
|
||||
The *property* setting can be *id* or *image* or *mol*. For *id*, the
|
||||
@ -90,7 +90,7 @@ keyword/value settings are given below.
|
||||
|
||||
----------
|
||||
|
||||
*Property id*
|
||||
Property: *id*
|
||||
|
||||
Reset atom IDs for the entire system, including all the global IDs
|
||||
stored for bond, angle, dihedral, improper topology data. This will
|
||||
@ -146,7 +146,7 @@ processor have consecutive IDs, as the :doc:`create_atoms
|
||||
|
||||
----------
|
||||
|
||||
*Property image*
|
||||
Property: *image*
|
||||
|
||||
Reset the image flags of atoms so that at least one atom in each
|
||||
molecule has an image flag of 0. Molecular topology is respected so
|
||||
@ -191,7 +191,7 @@ flags.
|
||||
|
||||
----------
|
||||
|
||||
*Property mol*
|
||||
Property: *mol*
|
||||
|
||||
Reset molecule IDs for a specified group of atoms based on current
|
||||
bond connectivity. This will typically create a new set of molecule
|
||||
@ -203,7 +203,7 @@ For purposes of this operation, molecules are identified by the current
|
||||
bond connectivity in the system, which may or may not be consistent with
|
||||
the current molecule IDs. A molecule in this context is a set of atoms
|
||||
connected to each other with explicit bonds. The specific algorithm
|
||||
used is the one of :doc:`compute fragment/atom <compute_cluster_atom>`
|
||||
used is the one of :doc:`compute fragment/atom <compute_cluster_atom>`.
|
||||
Once the molecules are identified and a new molecule ID computed for
|
||||
each, this command will update the current molecule ID for all atoms in
|
||||
the group with the new molecule ID. Note that if the group excludes
|
||||
@ -266,7 +266,7 @@ The *image* property can only be used when the atom style supports bonds.
|
||||
Related commands
|
||||
""""""""""""""""
|
||||
|
||||
:doc:`compute fragment/atom <compute_cluster_atom>`
|
||||
:doc:`compute fragment/atom <compute_cluster_atom>`,
|
||||
:doc:`fix bond/react <fix_bond_react>`,
|
||||
:doc:`fix bond/create <fix_bond_create>`,
|
||||
:doc:`fix bond/break <fix_bond_break>`,
|
||||
|
||||
@ -1,2 +1,3 @@
|
||||
*.csv
|
||||
*.txt
|
||||
*.lammpstrj
|
||||
|
||||
@ -17,14 +17,22 @@ q_ref = float(ref_line[3])
|
||||
inv11_ref = float(ref_line[4])
|
||||
inv12_ref = float(ref_line[5])
|
||||
b1_ref = float(ref_line[6])
|
||||
felec1_ref = float(ref_line[8])
|
||||
felyt1_ref = float(ref_line[10])
|
||||
press_ref = float(ref_line[12])
|
||||
|
||||
# out.csv
|
||||
with open(sys.argv[2]) as f:
|
||||
out_line = f.readlines()[-1].split(", ")
|
||||
e_out = float(out_line[0])
|
||||
q_out = float(out_line[1])
|
||||
press_out = float(out_line[2])
|
||||
|
||||
out_lines = [("energy", e_ref, e_out), ("charge", q_ref, q_out)]
|
||||
out_lines = [
|
||||
("energy", e_ref, e_out),
|
||||
("charge", q_ref, q_out),
|
||||
("pressure", press_ref, press_out),
|
||||
]
|
||||
|
||||
# vec.csv
|
||||
vec_file = "vec.csv"
|
||||
@ -44,6 +52,14 @@ if op.isfile(inv_file):
|
||||
inv12_out = float(inv_line[1])
|
||||
out_lines.append(("inv11", inv11_ref, inv11_out))
|
||||
|
||||
# forces.lammpstrj
|
||||
force_file = "forces.lammpstrj"
|
||||
with open(force_file) as f:
|
||||
lines = f.readlines()[9:]
|
||||
for name, i, f_ref in [("felec1", "1", felec1_ref), ("felyt1", "3", felyt1_ref)]:
|
||||
f_out = next(float(y[3]) for x in lines if (y := x.split()) and y[0] == i)
|
||||
out_lines.append((name, f_ref, f_out))
|
||||
|
||||
lines = []
|
||||
for label, ref, out in out_lines:
|
||||
error = rel_error(out, ref)
|
||||
|
||||
@ -8,7 +8,7 @@ thermo_style custom step pe c_qbot c_qtop
|
||||
fix feta all property/atom d_eta ghost on
|
||||
set group bot d_eta 0.5
|
||||
set group top d_eta 3.0
|
||||
fix conp bot electrode/conp 0 2 couple top 1 symm on eta d_eta algo cg 1e-6
|
||||
fix conp bot electrode/conp 0 NULL couple top 1 symm on eta d_eta algo cg 1e-6
|
||||
|
||||
run 0
|
||||
|
||||
|
||||
@ -8,7 +8,7 @@ thermo_style custom step pe c_qbot c_qtop
|
||||
fix feta all property/atom d_eta ghost on
|
||||
set group bot d_eta 0.5
|
||||
set group top d_eta 3.0
|
||||
fix conp bot electrode/conp 0 2 couple top 1 symm on eta d_eta write_inv inv.csv write_vec vec.csv
|
||||
fix conp bot electrode/conp 0 NULL couple top 1 symm on eta d_eta write_inv inv.csv write_vec vec.csv
|
||||
|
||||
run 0
|
||||
|
||||
|
||||
@ -1,12 +1,17 @@
|
||||
#!/usr/bin/env python3
|
||||
|
||||
import time
|
||||
|
||||
import numpy as np
|
||||
from scipy.special import erf
|
||||
|
||||
SQRT2 = np.sqrt(2)
|
||||
SQRTPI_INV = 1 / np.sqrt(np.pi)
|
||||
COULOMB = 332.06371 # Coulomb constant in Lammps 'real' units
|
||||
QE2F = 23.060549
|
||||
NKTV2P = 68568.415 # pressure in 'real' units
|
||||
LENGTH = 10000 # convergence parameter
|
||||
LZ = 20
|
||||
|
||||
|
||||
def lattice(length):
|
||||
@ -26,6 +31,25 @@ def b_element(r, q, eta):
|
||||
return q * erf(eta * r) / r
|
||||
|
||||
|
||||
def force_gauss(r, qq, eta):
|
||||
etar = eta * r
|
||||
return (qq / np.square(r)) * (
|
||||
erf(etar) - 2 * etar * SQRTPI_INV * np.exp(-np.square(etar))
|
||||
)
|
||||
|
||||
|
||||
def force_point(r, qq):
|
||||
return qq / np.square(r)
|
||||
|
||||
|
||||
def force_component(dx, d, qq, eta=None):
|
||||
if eta:
|
||||
return np.sum(dx / d * force_gauss(d, qq, eta))
|
||||
else:
|
||||
return np.sum(dx / d * force_point(d, qq))
|
||||
|
||||
|
||||
time_start = time.perf_counter()
|
||||
a = 1 # nearest neighbor distance i.e. lattice constant / sqrt(2)
|
||||
x_elec = [-2, 2]
|
||||
x_elyt = [-1, 1]
|
||||
@ -36,8 +60,20 @@ v = np.array([-0.5, 0.5]) * (QE2F / COULOMB)
|
||||
# distances to images within electrode and to opposite electrode
|
||||
distances = a * np.linalg.norm(lattice(LENGTH), axis=1)
|
||||
opposite_distances = np.sqrt(np.square(distances) + distance_plates**2)
|
||||
image_distances = []
|
||||
for x in x_elec:
|
||||
image_distances.append([])
|
||||
for y in x_elyt:
|
||||
image_distances[-1].append(np.sqrt(np.square(distances) + np.abs(y - x) ** 2))
|
||||
image_elyt_distances = [[None for _ in range(len(x_elyt))] for _ in range(len(x_elyt))]
|
||||
for i, (xi, qi) in enumerate(zip(x_elyt, q_elyt)):
|
||||
for j, (xj, qj) in list(enumerate(zip(x_elyt, q_elyt)))[i + 1 :]:
|
||||
image_elyt_distances[i][j] = np.sqrt(
|
||||
np.square(distances) + np.abs(xj - xi) ** 2
|
||||
)
|
||||
|
||||
for name, eta_elec in [("", [2.0, 2.0]), ("_eta_mix", [0.5, 3.0])]:
|
||||
# for name, eta_elec in [("", [2.0, 2.0])]:
|
||||
eta_mix = np.prod(eta_elec) / np.sqrt(np.sum(np.square(eta_elec)))
|
||||
# self interaction and within original box
|
||||
A_11 = np.sqrt(2 / np.pi) * eta_elec[0]
|
||||
@ -55,22 +91,18 @@ for name, eta_elec in [("", [2.0, 2.0]), ("_eta_mix", [0.5, 3.0])]:
|
||||
|
||||
# electrode-electrolyte interaction
|
||||
b = []
|
||||
for x, eta in zip(x_elec, eta_elec):
|
||||
for i, (x, eta) in enumerate(zip(x_elec, eta_elec)):
|
||||
bi = 0
|
||||
for y, q in zip(x_elyt, q_elyt):
|
||||
d = abs(y - x)
|
||||
bi += b_element(d, q, eta)
|
||||
image_distances = np.sqrt(np.square(distances) + d**2)
|
||||
bi += 4 * np.sum(b_element(image_distances, q, eta))
|
||||
for j, (y, q) in enumerate(zip(x_elyt, q_elyt)):
|
||||
bi += b_element(np.abs(y - x), q, eta)
|
||||
bi += 4 * np.sum(b_element(image_distances[i][j], q, eta))
|
||||
b.append(bi)
|
||||
b = np.array(b)
|
||||
|
||||
# electrolyte-electrolyte energy
|
||||
elyt_11 = 4 * np.sum(1 / distances)
|
||||
distance_elyt = x_elyt[1] - x_elyt[0]
|
||||
elyt_12 = 1 / distance_elyt + 4 * np.sum(
|
||||
1 / np.sqrt(np.square(distances) + distance_elyt**2)
|
||||
)
|
||||
elyt_12 = 1 / distance_elyt + 4 * np.sum(1 / image_elyt_distances[0][1])
|
||||
elyt = np.array([[elyt_11, elyt_12], [elyt_12, elyt_11]])
|
||||
energy_elyt = 0.5 * np.dot(q_elyt, np.dot(elyt, q_elyt))
|
||||
|
||||
@ -78,9 +110,48 @@ for name, eta_elec in [("", [2.0, 2.0]), ("_eta_mix", [0.5, 3.0])]:
|
||||
q = np.dot(inv, v - b)
|
||||
energy = COULOMB * (0.5 * np.dot(q, np.dot(A, q)) + np.dot(b, q) + energy_elyt)
|
||||
|
||||
# forces in out-of-plane direction
|
||||
f_elec = np.zeros(len(x_elec))
|
||||
f_elyt = np.zeros(len(x_elyt))
|
||||
# electrode-electrode
|
||||
dx = x_elec[1] - x_elec[0]
|
||||
fij_box = force_component(dx, np.abs(dx), q[0] * q[1], eta_mix)
|
||||
fij_img = 4 * force_component(dx, opposite_distances, q[0] * q[1], eta_mix)
|
||||
f_elec[0] -= fij_box + fij_img
|
||||
f_elec[1] += fij_box + fij_img
|
||||
# electrode-electrolyte
|
||||
for i, (xi, qi, etai) in enumerate(zip(x_elec, q, eta_elec)):
|
||||
for j, (xj, qj) in enumerate(zip(x_elyt, q_elyt)):
|
||||
dx = xj - xi
|
||||
fij_box = force_component(dx, np.abs(dx), qi * qj, etai)
|
||||
fij_img = 4 * force_component(dx, image_distances[i][j], qi * qj, etai)
|
||||
f_elec[i] -= fij_box + fij_img
|
||||
f_elyt[j] += fij_box + fij_img
|
||||
# electrolyte-electrolyte
|
||||
for i, (xi, qi) in enumerate(zip(x_elyt, q_elyt)):
|
||||
for j, (xj, qj) in list(enumerate(zip(x_elyt, q_elyt)))[i + 1 :]:
|
||||
dx = xj - xi
|
||||
fij_box = force_component(dx, np.abs(dx), qi * qj)
|
||||
fij_img = 4 * force_component(dx, image_elyt_distances[i][j], qi * qj)
|
||||
f_elyt[i] -= fij_img + fij_box
|
||||
f_elyt[j] += fij_img + fij_box
|
||||
# force units
|
||||
assert np.abs(np.sum(f_elec) + np.sum(f_elyt)) < 1e-8
|
||||
f_elec *= COULOMB
|
||||
f_elyt *= COULOMB
|
||||
|
||||
# Virial
|
||||
volume = a**2 * LZ
|
||||
virial = 0.0
|
||||
for x, f in [(x_elec, f_elec), (x_elyt, f_elyt)]:
|
||||
virial += np.dot(x, f)
|
||||
pressure = NKTV2P * virial / volume
|
||||
|
||||
with open(f"plate_cap{name}.csv", "w") as f:
|
||||
f.write(
|
||||
"length, energy / kcal/mol, q1 / e, q2 / e, inv11 / A, inv12 / A, b1 / e/A, b2 / e/A\n"
|
||||
"length, energy / kcal/mol, q1 / e, q2 / e, inv11 / A, inv12 / A"
|
||||
+ ", b1 / e/A, b2 / e/A, felec1 / kcal/mol/A, felec2 / kcal/mol/A"
|
||||
+ ", felyt1 / kcal/mol/A, felyt2 / kcal/mol/A, press\n"
|
||||
)
|
||||
f.write(
|
||||
", ".join(
|
||||
@ -93,7 +164,14 @@ for name, eta_elec in [("", [2.0, 2.0]), ("_eta_mix", [0.5, 3.0])]:
|
||||
f"{inv[0, 1]:.10f}",
|
||||
f"{b[0]:.8f}",
|
||||
f"{b[1]:.8f}",
|
||||
f"{f_elec[0]:.5f}",
|
||||
f"{f_elec[1]:.5f}",
|
||||
f"{f_elyt[0]:.5f}",
|
||||
f"{f_elyt[1]:.5f}",
|
||||
f"{pressure:.2f}",
|
||||
]
|
||||
)
|
||||
+ "\n"
|
||||
)
|
||||
time_end = time.perf_counter()
|
||||
print(f"{time_end - time_start:0.4f} seconds")
|
||||
|
||||
@ -19,4 +19,8 @@ compute qtop top reduce sum v_q
|
||||
compute compute_pe all pe
|
||||
variable vpe equal c_compute_pe
|
||||
variable charge equal c_qtop
|
||||
fix fxprint all print 1 "${vpe}, ${charge}" file "out.csv"
|
||||
compute press all pressure NULL virial
|
||||
variable p3 equal c_press[3]
|
||||
fix fxprint all print 1 "${vpe}, ${charge}, ${p3}" file "out.csv"
|
||||
|
||||
dump dump_forces all custom 1 forces.lammpstrj id fx fy fz
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
@ -42,7 +42,11 @@ fix fxforce_au gold setforce 0.0 0.0 0.0
|
||||
|
||||
# equilibrate z-coordinate of upper electrode while keeping the electrode rigid
|
||||
fix fxforce_wa wall setforce 0.0 0.0 NULL
|
||||
fix fxpressure wall aveforce 0 0 -0.005246 # atomspheric pressure: area/force->nktv2p
|
||||
variable atm equal 1/68568.415 # 1/force->nktv2p
|
||||
variable area equal (xhi-xlo)*(yhi-ylo)
|
||||
variable wall_force equal -v_atm*v_area/count(wall)
|
||||
print "Wall force per atom: ${wall_force}"
|
||||
fix fxpressure wall aveforce 0 0 ${wall_force} # atomspheric pressure: area/force->nktv2p
|
||||
fix fxdrag wall viscous 100
|
||||
fix fxrigid wall rigid/nve single
|
||||
|
||||
|
||||
@ -1,4 +1,5 @@
|
||||
LAMMPS (3 Nov 2022)
|
||||
LAMMPS (7 Feb 2024 - Development - patch_7Feb2024_update1-217-g1909233c69-modified)
|
||||
using 1 OpenMP thread(s) per MPI task
|
||||
# The intention is to find the average position of one wall at atmospheric
|
||||
# pressure. The output is the wall position over time which can be used to
|
||||
# find the average position for a run with fixed wall position.
|
||||
@ -40,8 +41,8 @@ Finding 1-2 1-3 1-4 neighbors ...
|
||||
1 = max # of 1-3 neighbors
|
||||
1 = max # of 1-4 neighbors
|
||||
2 = max # of special neighbors
|
||||
special bonds CPU = 0.001 seconds
|
||||
read_data CPU = 0.011 seconds
|
||||
special bonds CPU = 0.000 seconds
|
||||
read_data CPU = 0.012 seconds
|
||||
|
||||
# ----------------- Settings Section -----------------
|
||||
|
||||
@ -77,7 +78,13 @@ fix fxforce_au gold setforce 0.0 0.0 0.0
|
||||
|
||||
# equilibrate z-coordinate of upper electrode while keeping the electrode rigid
|
||||
fix fxforce_wa wall setforce 0.0 0.0 NULL
|
||||
fix fxpressure wall aveforce 0 0 -0.005246 # atomspheric pressure: area/force->nktv2p
|
||||
variable atm equal 1/68568.415 # 1/force->nktv2p
|
||||
variable area equal (xhi-xlo)*(yhi-ylo)
|
||||
variable wall_force equal -v_atm*v_area/count(wall)
|
||||
print "Wall force per atom: ${wall_force}"
|
||||
Wall force per atom: -0.000109285996244287
|
||||
fix fxpressure wall aveforce 0 0 ${wall_force} # atomspheric pressure: area/force->nktv2p
|
||||
fix fxpressure wall aveforce 0 0 -0.000109285996244287
|
||||
fix fxdrag wall viscous 100
|
||||
fix fxrigid wall rigid/nve single
|
||||
1 rigid bodies with 48 atoms
|
||||
@ -134,7 +141,7 @@ PPPM/electrode initialization ...
|
||||
stencil order = 5
|
||||
estimated absolute RMS force accuracy = 0.02930901
|
||||
estimated relative force accuracy = 8.8263214e-05
|
||||
using double precision MKL FFT
|
||||
using double precision FFTW3
|
||||
3d grid and FFT values/proc = 15884 6480
|
||||
Generated 6 of 6 mixed pair_coeff terms from arithmetic mixing rule
|
||||
Neighbor list info ...
|
||||
@ -157,54 +164,54 @@ Neighbor list info ...
|
||||
Per MPI rank memory allocation (min/avg/max) = 11.7 | 11.7 | 11.7 Mbytes
|
||||
Step c_temp_mobile c_qwa c_qau v_top_wall
|
||||
0 303.38967 -0.042963484 0.042963484 21.4018
|
||||
5000 285.08828 -0.26105255 0.26105255 25.155629
|
||||
10000 323.19176 -0.26264003 0.26264003 24.541676
|
||||
15000 310.479 -0.27318148 0.27318148 23.141522
|
||||
20000 295.18544 -0.11313444 0.11313444 23.828735
|
||||
25000 295.38607 -0.25433086 0.25433086 23.673314
|
||||
30000 288.0613 -0.30099901 0.30099901 23.438086
|
||||
35000 278.5591 -0.15823576 0.15823576 24.311915
|
||||
40000 303.95751 -0.19941381 0.19941381 23.69594
|
||||
45000 279.026 -0.1659962 0.1659962 23.588604
|
||||
50000 298.79278 -0.28866703 0.28866703 23.372508
|
||||
55000 301.03353 -0.078370381 0.078370381 23.192985
|
||||
60000 306.77965 -0.12807205 0.12807205 23.968574
|
||||
65000 309.86008 -0.27162663 0.27162663 23.616704
|
||||
70000 287.31116 -0.029751882 0.029751882 23.667495
|
||||
75000 312.48654 -0.10759866 0.10759866 23.504105
|
||||
80000 309.94267 -0.2558548 0.2558548 23.810576
|
||||
85000 328.04389 -0.1575471 0.1575471 24.013437
|
||||
90000 302.9806 -0.032002164 0.032002164 24.264432
|
||||
95000 294.20804 -0.27797238 0.27797238 23.291758
|
||||
100000 307.63019 -0.19047448 0.19047448 23.632147
|
||||
Loop time of 530.844 on 1 procs for 100000 steps with 726 atoms
|
||||
5000 311.85363 0.03543775 -0.03543775 24.79665
|
||||
10000 285.91321 -0.16873703 0.16873703 23.103088
|
||||
15000 295.39476 -0.44424612 0.44424612 23.767107
|
||||
20000 296.12969 -0.14120993 0.14120993 23.96361
|
||||
25000 306.59629 -0.29333182 0.29333182 23.884488
|
||||
30000 297.98559 -0.10749684 0.10749684 23.73316
|
||||
35000 297.98503 -0.11809975 0.11809975 23.984669
|
||||
40000 300.26292 -0.32784184 0.32784184 23.462748
|
||||
45000 295.68441 -0.25940165 0.25940165 23.516403
|
||||
50000 315.12883 -0.36037614 0.36037614 23.627879
|
||||
55000 290.55151 -0.0032838106 0.0032838106 23.684931
|
||||
60000 316.4625 -0.17245368 0.17245368 24.126883
|
||||
65000 296.79343 -0.054061851 0.054061851 23.695094
|
||||
70000 305.99923 -0.11363801 0.11363801 23.55476
|
||||
75000 297.40131 -0.27054153 0.27054153 23.928994
|
||||
80000 306.54811 -0.25409719 0.25409719 23.869448
|
||||
85000 303.95231 -0.17895561 0.17895561 23.658833
|
||||
90000 313.43739 -0.059036514 0.059036514 23.36056
|
||||
95000 290.3077 -0.31394478 0.31394478 23.885538
|
||||
100000 297.5156 -0.30730083 0.30730083 23.511674
|
||||
Loop time of 1586.06 on 1 procs for 100000 steps with 726 atoms
|
||||
|
||||
Performance: 32.552 ns/day, 0.737 hours/ns, 188.379 timesteps/s, 136.763 katom-step/s
|
||||
100.0% CPU use with 1 MPI tasks x no OpenMP threads
|
||||
Performance: 10.895 ns/day, 2.203 hours/ns, 63.049 timesteps/s, 45.774 katom-step/s
|
||||
99.6% CPU use with 1 MPI tasks x 1 OpenMP threads
|
||||
|
||||
MPI task timing breakdown:
|
||||
Section | min time | avg time | max time |%varavg| %total
|
||||
---------------------------------------------------------------
|
||||
Pair | 190.47 | 190.47 | 190.47 | 0.0 | 35.88
|
||||
Bond | 0.10754 | 0.10754 | 0.10754 | 0.0 | 0.02
|
||||
Kspace | 73.179 | 73.179 | 73.179 | 0.0 | 13.79
|
||||
Neigh | 24.209 | 24.209 | 24.209 | 0.0 | 4.56
|
||||
Comm | 1.6857 | 1.6857 | 1.6857 | 0.0 | 0.32
|
||||
Output | 0.0016861 | 0.0016861 | 0.0016861 | 0.0 | 0.00
|
||||
Modify | 240.23 | 240.23 | 240.23 | 0.0 | 45.26
|
||||
Other | | 0.9595 | | | 0.18
|
||||
Pair | 460.91 | 460.91 | 460.91 | 0.0 | 29.06
|
||||
Bond | 0.047873 | 0.047873 | 0.047873 | 0.0 | 0.00
|
||||
Kspace | 341.4 | 341.4 | 341.4 | 0.0 | 21.53
|
||||
Neigh | 52.868 | 52.868 | 52.868 | 0.0 | 3.33
|
||||
Comm | 5.2321 | 5.2321 | 5.2321 | 0.0 | 0.33
|
||||
Output | 0.00099102 | 0.00099102 | 0.00099102 | 0.0 | 0.00
|
||||
Modify | 724.63 | 724.63 | 724.63 | 0.0 | 45.69
|
||||
Other | | 0.9741 | | | 0.06
|
||||
|
||||
Nlocal: 726 ave 726 max 726 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
Nghost: 2335 ave 2335 max 2335 min
|
||||
Nghost: 2336 ave 2336 max 2336 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
Neighs: 120271 ave 120271 max 120271 min
|
||||
Neighs: 120321 ave 120321 max 120321 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
|
||||
Total # of neighbors = 120271
|
||||
Ave neighs/atom = 165.66253
|
||||
Total # of neighbors = 120321
|
||||
Ave neighs/atom = 165.7314
|
||||
Ave special neighs/atom = 1.7355372
|
||||
Neighbor list builds = 7722
|
||||
Neighbor list builds = 7670
|
||||
Dangerous builds = 0
|
||||
write_data "data.piston.final"
|
||||
System init for write_data ...
|
||||
@ -213,11 +220,11 @@ PPPM/electrode initialization ...
|
||||
G vector (1/distance) = 0.32814871
|
||||
grid = 12 15 36
|
||||
stencil order = 5
|
||||
estimated absolute RMS force accuracy = 0.029311365
|
||||
estimated relative force accuracy = 8.8270304e-05
|
||||
using double precision MKL FFT
|
||||
estimated absolute RMS force accuracy = 0.029311329
|
||||
estimated relative force accuracy = 8.8270197e-05
|
||||
using double precision FFTW3
|
||||
3d grid and FFT values/proc = 15884 6480
|
||||
Generated 6 of 6 mixed pair_coeff terms from arithmetic mixing rule
|
||||
|
||||
Average conjugate gradient steps: 1.981
|
||||
Total wall time: 0:08:50
|
||||
Total wall time: 0:26:26
|
||||
@ -1,4 +1,5 @@
|
||||
LAMMPS (3 Nov 2022)
|
||||
LAMMPS (7 Feb 2024 - Development - patch_7Feb2024_update1-217-g1909233c69-modified)
|
||||
using 1 OpenMP thread(s) per MPI task
|
||||
# The intention is to find the average position of one wall at atmospheric
|
||||
# pressure. The output is the wall position over time which can be used to
|
||||
# find the average position for a run with fixed wall position.
|
||||
@ -41,8 +42,8 @@ Finding 1-2 1-3 1-4 neighbors ...
|
||||
1 = max # of 1-3 neighbors
|
||||
1 = max # of 1-4 neighbors
|
||||
2 = max # of special neighbors
|
||||
special bonds CPU = 0.001 seconds
|
||||
read_data CPU = 0.017 seconds
|
||||
special bonds CPU = 0.000 seconds
|
||||
read_data CPU = 0.012 seconds
|
||||
|
||||
# ----------------- Settings Section -----------------
|
||||
|
||||
@ -66,7 +67,7 @@ Finding SHAKE clusters ...
|
||||
0 = # of size 3 clusters
|
||||
0 = # of size 4 clusters
|
||||
210 = # of frozen angles
|
||||
find clusters CPU = 0.002 seconds
|
||||
find clusters CPU = 0.000 seconds
|
||||
pair_modify mix arithmetic
|
||||
|
||||
# ----------------- Run Section -----------------
|
||||
@ -78,7 +79,13 @@ fix fxforce_au gold setforce 0.0 0.0 0.0
|
||||
|
||||
# equilibrate z-coordinate of upper electrode while keeping the electrode rigid
|
||||
fix fxforce_wa wall setforce 0.0 0.0 NULL
|
||||
fix fxpressure wall aveforce 0 0 -0.005246 # atomspheric pressure: area/force->nktv2p
|
||||
variable atm equal 1/68568.415 # 1/force->nktv2p
|
||||
variable area equal (xhi-xlo)*(yhi-ylo)
|
||||
variable wall_force equal -v_atm*v_area/count(wall)
|
||||
print "Wall force per atom: ${wall_force}"
|
||||
Wall force per atom: -0.000109285996244287
|
||||
fix fxpressure wall aveforce 0 0 ${wall_force} # atomspheric pressure: area/force->nktv2p
|
||||
fix fxpressure wall aveforce 0 0 -0.000109285996244287
|
||||
fix fxdrag wall viscous 100
|
||||
fix fxrigid wall rigid/nve single
|
||||
1 rigid bodies with 48 atoms
|
||||
@ -135,7 +142,7 @@ PPPM/electrode initialization ...
|
||||
stencil order = 5
|
||||
estimated absolute RMS force accuracy = 0.02930901
|
||||
estimated relative force accuracy = 8.8263214e-05
|
||||
using double precision MKL FFT
|
||||
using double precision FFTW3
|
||||
3d grid and FFT values/proc = 8512 2880
|
||||
Generated 6 of 6 mixed pair_coeff terms from arithmetic mixing rule
|
||||
Neighbor list info ...
|
||||
@ -158,54 +165,54 @@ Neighbor list info ...
|
||||
Per MPI rank memory allocation (min/avg/max) = 10.06 | 10.22 | 10.41 Mbytes
|
||||
Step c_temp_mobile c_qwa c_qau v_top_wall
|
||||
0 303.38967 -0.042963484 0.042963484 21.4018
|
||||
5000 292.03027 -0.19040435 0.19040435 24.581338
|
||||
10000 309.52764 -0.48308301 0.48308301 23.776985
|
||||
15000 295.00243 -0.16591109 0.16591109 23.672038
|
||||
20000 293.5536 -0.086669084 0.086669084 23.426455
|
||||
25000 303.0079 -0.16488112 0.16488112 23.862966
|
||||
30000 306.31463 -0.23192653 0.23192653 23.819882
|
||||
35000 303.66268 -0.2317907 0.2317907 23.495344
|
||||
40000 301.39435 -0.34661329 0.34661329 23.657835
|
||||
45000 291.61205 -0.30539427 0.30539427 23.437303
|
||||
50000 298.65319 -0.096107034 0.096107034 23.57809
|
||||
55000 282.65069 -0.14943539 0.14943539 23.823728
|
||||
60000 310.64182 -0.17418813 0.17418813 23.286959
|
||||
65000 308.47141 -0.02075662 0.02075662 23.91313
|
||||
70000 292.5186 -0.080163162 0.080163162 23.96283
|
||||
75000 270.13928 -0.029528648 0.029528648 23.488972
|
||||
80000 322.10914 0.030761045 -0.030761045 23.47592
|
||||
85000 310.60347 -0.24069996 0.24069996 23.987091
|
||||
90000 294.35695 -0.070458235 0.070458235 23.397929
|
||||
95000 308.69043 -0.2652581 0.2652581 23.473813
|
||||
100000 318.71883 0.024035956 -0.024035956 23.449863
|
||||
Loop time of 590.232 on 4 procs for 100000 steps with 726 atoms
|
||||
5000 291.6303 -0.1820085 0.1820085 24.641399
|
||||
10000 299.42886 -0.19823095 0.19823095 23.820522
|
||||
15000 288.23071 -0.065261869 0.065261869 23.360845
|
||||
20000 299.4644 -0.042993777 0.042993777 23.987554
|
||||
25000 304.26497 -0.15665293 0.15665293 23.729006
|
||||
30000 292.29674 -0.25142779 0.25142779 23.960725
|
||||
35000 295.57492 -0.01269228 0.01269228 23.445383
|
||||
40000 303.38438 -0.13941727 0.13941727 23.517483
|
||||
45000 302.211 -0.19589892 0.19589892 23.704043
|
||||
50000 281.64939 -0.18057298 0.18057298 23.542137
|
||||
55000 274.90565 -0.15453379 0.15453379 23.734347
|
||||
60000 290.70459 -0.27977436 0.27977436 23.835365
|
||||
65000 293.42241 -0.2454241 0.2454241 23.59269
|
||||
70000 295.20229 -0.041314995 0.041314995 23.73856
|
||||
75000 297.79519 -0.11231755 0.11231755 23.57262
|
||||
80000 285.17858 -0.070796508 0.070796508 23.817135
|
||||
85000 311.71609 -0.068920177 0.068920177 23.861127
|
||||
90000 287.80446 -0.19183387 0.19183387 23.369393
|
||||
95000 309.43345 -0.15238671 0.15238671 23.597792
|
||||
100000 294.12422 -0.14284353 0.14284353 23.526286
|
||||
Loop time of 876.546 on 4 procs for 100000 steps with 726 atoms
|
||||
|
||||
Performance: 29.277 ns/day, 0.820 hours/ns, 169.425 timesteps/s, 123.003 katom-step/s
|
||||
72.5% CPU use with 4 MPI tasks x no OpenMP threads
|
||||
Performance: 19.714 ns/day, 1.217 hours/ns, 114.084 timesteps/s, 82.825 katom-step/s
|
||||
98.6% CPU use with 4 MPI tasks x 1 OpenMP threads
|
||||
|
||||
MPI task timing breakdown:
|
||||
Section | min time | avg time | max time |%varavg| %total
|
||||
---------------------------------------------------------------
|
||||
Pair | 57.391 | 75.867 | 96.292 | 212.1 | 12.85
|
||||
Bond | 0.10177 | 0.11042 | 0.12415 | 2.7 | 0.02
|
||||
Kspace | 102.79 | 123.16 | 141.5 | 165.7 | 20.87
|
||||
Neigh | 12.808 | 12.895 | 12.982 | 2.3 | 2.18
|
||||
Comm | 18.885 | 19.973 | 21.064 | 24.0 | 3.38
|
||||
Output | 0.0022573 | 0.0022749 | 0.0023225 | 0.1 | 0.00
|
||||
Modify | 355.89 | 356.74 | 357.61 | 4.2 | 60.44
|
||||
Other | | 1.478 | | | 0.25
|
||||
Pair | 123.63 | 171.23 | 215.73 | 336.6 | 19.53
|
||||
Bond | 0.068261 | 0.075883 | 0.081822 | 1.9 | 0.01
|
||||
Kspace | 187.59 | 231.71 | 279.01 | 287.1 | 26.43
|
||||
Neigh | 29.28 | 29.462 | 29.637 | 2.5 | 3.36
|
||||
Comm | 12.544 | 13.731 | 14.929 | 29.1 | 1.57
|
||||
Output | 0.0010182 | 0.0014585 | 0.0016071 | 0.7 | 0.00
|
||||
Modify | 428.74 | 429.25 | 429.74 | 2.3 | 48.97
|
||||
Other | | 1.092 | | | 0.12
|
||||
|
||||
Nlocal: 181.5 ave 207 max 169 min
|
||||
Histogram: 2 0 1 0 0 0 0 0 0 1
|
||||
Nghost: 1961.5 ave 1984 max 1926 min
|
||||
Histogram: 1 0 0 0 0 0 1 0 1 1
|
||||
Neighs: 30051 ave 41646 max 20775 min
|
||||
Histogram: 1 1 0 0 0 0 1 0 0 1
|
||||
Nlocal: 181.5 ave 195 max 166 min
|
||||
Histogram: 1 1 0 0 0 0 0 0 0 2
|
||||
Nghost: 1955.5 ave 1978 max 1931 min
|
||||
Histogram: 1 0 0 0 1 0 1 0 0 1
|
||||
Neighs: 30343 ave 39847 max 20428 min
|
||||
Histogram: 2 0 0 0 0 0 0 0 0 2
|
||||
|
||||
Total # of neighbors = 120204
|
||||
Ave neighs/atom = 165.57025
|
||||
Total # of neighbors = 121372
|
||||
Ave neighs/atom = 167.17906
|
||||
Ave special neighs/atom = 1.7355372
|
||||
Neighbor list builds = 7663
|
||||
Neighbor list builds = 7698
|
||||
Dangerous builds = 0
|
||||
write_data "data.piston.final"
|
||||
System init for write_data ...
|
||||
@ -214,11 +221,11 @@ PPPM/electrode initialization ...
|
||||
G vector (1/distance) = 0.32814871
|
||||
grid = 12 15 36
|
||||
stencil order = 5
|
||||
estimated absolute RMS force accuracy = 0.029311028
|
||||
estimated relative force accuracy = 8.8269289e-05
|
||||
using double precision MKL FFT
|
||||
estimated absolute RMS force accuracy = 0.029310954
|
||||
estimated relative force accuracy = 8.8269069e-05
|
||||
using double precision FFTW3
|
||||
3d grid and FFT values/proc = 8512 2880
|
||||
Generated 6 of 6 mixed pair_coeff terms from arithmetic mixing rule
|
||||
|
||||
Average conjugate gradient steps: 1.982
|
||||
Total wall time: 0:09:50
|
||||
Average conjugate gradient steps: 1.981
|
||||
Total wall time: 0:14:36
|
||||
@ -79,6 +79,9 @@ FixElectrodeConp::FixElectrodeConp(LAMMPS *lmp, int narg, char **arg) :
|
||||
potential_i(nullptr), potential_iele(nullptr)
|
||||
{
|
||||
if (lmp->citeme) lmp->citeme->add(cite_fix_electrode);
|
||||
if (atom->map_style == Atom::MAP_NONE)
|
||||
error->all(FLERR, "Fix {} requires an atom map, see atom_modify", style);
|
||||
|
||||
// fix.h output flags
|
||||
scalar_flag = 1;
|
||||
vector_flag = 1;
|
||||
@ -86,6 +89,9 @@ FixElectrodeConp::FixElectrodeConp(LAMMPS *lmp, int narg, char **arg) :
|
||||
extvector = 0;
|
||||
extarray = 0;
|
||||
|
||||
virial_global_flag = 1; // use virials of this fix
|
||||
thermo_virial = 1; // set vflags for v_tally
|
||||
|
||||
bool default_algo = true;
|
||||
algo = Algo::MATRIX_INV;
|
||||
matrix_algo = true;
|
||||
@ -123,7 +129,8 @@ FixElectrodeConp::FixElectrodeConp(LAMMPS *lmp, int narg, char **arg) :
|
||||
} else
|
||||
group_psi_const[0] = utils::numeric(FLERR, arg[3], false, lmp);
|
||||
char *eta_str = arg[4];
|
||||
eta = utils::numeric(FLERR, eta_str, false, lmp);
|
||||
bool etanull = (strcmp(eta_str, "NULL") == 0);
|
||||
if (!etanull) eta = utils::numeric(FLERR, eta_str, false, lmp);
|
||||
int iarg = 5;
|
||||
while (iarg < narg) {
|
||||
if ((strcmp(arg[iarg], "couple") == 0)) {
|
||||
@ -214,7 +221,7 @@ FixElectrodeConp::FixElectrodeConp(LAMMPS *lmp, int narg, char **arg) :
|
||||
qtotal_var_style = VarStyle::CONST;
|
||||
}
|
||||
} else if ((strcmp(arg[iarg], "eta") == 0)) {
|
||||
if (iarg + 2 > narg) error->all(FLERR, "Need two arguments after eta command");
|
||||
if (iarg + 2 > narg) error->all(FLERR, "Need one argument after eta command");
|
||||
etaflag = true;
|
||||
int is_double, cols, ghost;
|
||||
eta_index = atom->find_custom_ghost(arg[++iarg] + 2, is_double, cols, ghost);
|
||||
@ -240,6 +247,7 @@ FixElectrodeConp::FixElectrodeConp(LAMMPS *lmp, int narg, char **arg) :
|
||||
if (qtotal_var_style != VarStyle::UNSET) {
|
||||
if (symm) error->all(FLERR, "{} cannot use qtotal keyword with symm on", this->style);
|
||||
}
|
||||
if (etanull && !etaflag) error->all(FLERR, "If eta is NULL the eta keyword must be used");
|
||||
|
||||
// computatonal potential
|
||||
group_psi = std::vector<double>(groups.size());
|
||||
@ -346,7 +354,7 @@ int FixElectrodeConp::modify_param(int narg, char **arg)
|
||||
MPI_Allreduce(MPI_IN_PLACE, &in_ele, 1, MPI_INT, MPI_SUM, world);
|
||||
if (in_ele == 0) error->all(FLERR, "No atoms of type in electrode");
|
||||
MPI_Allreduce(MPI_IN_PLACE, ¬_in_ele, 1, MPI_INT, MPI_SUM, world);
|
||||
if (not_in_ele)
|
||||
if (not_in_ele && (comm->me == 0))
|
||||
error->warning(FLERR,
|
||||
"Not all atoms of type in electrode; Thomas-Fermi parameters will be ignored "
|
||||
"for electrolyte");
|
||||
@ -427,12 +435,11 @@ void FixElectrodeConp::init()
|
||||
}
|
||||
if (comm->me == 0)
|
||||
for (char *fix_id : integrate_ids)
|
||||
error->warning(
|
||||
FLERR,
|
||||
"Electrode atoms are integrated by fix {}, but fix electrode is using a matrix method. "
|
||||
"For "
|
||||
"mobile electrodes use the conjugate gradient algorithm without matrix ('algo cg').",
|
||||
fix_id);
|
||||
error->warning(FLERR,
|
||||
"Electrode atoms are integrated by fix {}, but fix electrode is using a "
|
||||
"matrix method. For mobile electrodes use the conjugate gradient algorithm "
|
||||
"without matrix ('algo cg').",
|
||||
fix_id);
|
||||
}
|
||||
|
||||
// check for package intel
|
||||
@ -478,7 +485,7 @@ void FixElectrodeConp::post_constructor()
|
||||
input->variable->set(fmt::format("{} equal f_{}[{}]", var_vtop, fixname, 1 + top_group));
|
||||
input->variable->set(fmt::format("{} equal (v_{}-v_{})/lz", var_efield, var_vbot, var_vtop));
|
||||
// check for other efields and warn if found
|
||||
if (modify->get_fix_by_style("^efield").size() > 0 && comm->me == 0)
|
||||
if ((modify->get_fix_by_style("^efield").size() > 0) && (comm->me == 0))
|
||||
error->warning(FLERR, "Other efield fixes found -- please make sure this is intended!");
|
||||
// call fix command:
|
||||
// fix [varstem]_efield all efield 0.0 0.0 [var_vdiff]/lz
|
||||
@ -638,13 +645,14 @@ void FixElectrodeConp::setup_post_neighbor()
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixElectrodeConp::setup_pre_reverse(int eflag, int /*vflag*/)
|
||||
void FixElectrodeConp::setup_pre_reverse(int eflag, int vflag)
|
||||
{
|
||||
if (pair->did_tally_callback() && (comm->me == 0))
|
||||
error->warning(FLERR, "Computation of virials in fix {} is incompatible with TALLY package", style);
|
||||
// correct forces for initial timestep
|
||||
gausscorr(eflag, true);
|
||||
ev_init(eflag, vflag);
|
||||
gausscorr(eflag, vflag, true);
|
||||
self_energy(eflag);
|
||||
// potential_energy(eflag); // not always part of the energy, depending on ensemble, therefore
|
||||
// removed
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
@ -674,7 +682,7 @@ void FixElectrodeConp::invert()
|
||||
void FixElectrodeConp::symmetrize()
|
||||
{
|
||||
// S matrix to enforce charge neutrality constraint
|
||||
if (read_inv && comm->me == 0)
|
||||
if (read_inv && (comm->me == 0))
|
||||
error->warning(FLERR,
|
||||
"Symmetrizing matrix from file. Make sure the provided matrix has not been "
|
||||
"symmetrized yet.");
|
||||
@ -755,11 +763,11 @@ void FixElectrodeConp::setup_pre_exchange() // create_taglist
|
||||
// if memory_usage > 0.5 GiB, warn with expected usage
|
||||
double mem_needed = memory_usage();
|
||||
mem_needed /= (1024 * 1024 * 1024); // convert to GiB
|
||||
if (mem_needed > 0.5 && comm->me == 0)
|
||||
if ((mem_needed > 0.5) && (comm->me == 0))
|
||||
error->warning(FLERR,
|
||||
"Please ensure there is sufficient memory for fix electrode "
|
||||
"Please ensure there is sufficient memory for fix {} "
|
||||
"(anticipated usage is at least {:.1f} GiB per proc)",
|
||||
mem_needed);
|
||||
style, mem_needed);
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
@ -771,12 +779,11 @@ void FixElectrodeConp::pre_force(int)
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixElectrodeConp::pre_reverse(int eflag, int /*vflag*/)
|
||||
void FixElectrodeConp::pre_reverse(int eflag, int vflag)
|
||||
{
|
||||
gausscorr(eflag, true);
|
||||
ev_init(eflag, vflag);
|
||||
gausscorr(eflag, vflag, true);
|
||||
self_energy(eflag);
|
||||
//potential_energy(eflag); // not always part of the energy, depending on ensemble, therefore
|
||||
// removed
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
@ -924,7 +931,7 @@ void FixElectrodeConp::update_charges()
|
||||
dot_old = dot_new;
|
||||
}
|
||||
recompute_potential(std::move(b), q_local);
|
||||
if (delta > cg_threshold && comm->me == 0) error->warning(FLERR, "CG threshold not reached");
|
||||
if ((delta > cg_threshold) && (comm->me == 0)) error->warning(FLERR, "CG threshold not reached");
|
||||
} else {
|
||||
error->all(FLERR, "This algorithm is not implemented, yet");
|
||||
}
|
||||
@ -1224,11 +1231,10 @@ double FixElectrodeConp::self_energy(int eflag)
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
double FixElectrodeConp::gausscorr(int eflag, bool fflag)
|
||||
double FixElectrodeConp::gausscorr(int eflag, int vflag, bool fflag)
|
||||
{
|
||||
// correction to short range interaction due to eta
|
||||
|
||||
int evflag = pair->evflag;
|
||||
double const qqrd2e = force->qqrd2e;
|
||||
int const nlocal = atom->nlocal;
|
||||
int *mask = atom->mask;
|
||||
@ -1294,13 +1300,11 @@ double FixElectrodeConp::gausscorr(int eflag, bool fflag)
|
||||
f[j][2] -= delz * fpair;
|
||||
}
|
||||
}
|
||||
|
||||
double ecoul = 0.;
|
||||
if (eflag) ecoul = -prefactor * erfc_etar;
|
||||
|
||||
if (evflag) {
|
||||
force->pair->ev_tally(i, j, nlocal, newton_pair, 0., ecoul, fpair, delx, dely, delz);
|
||||
if (eflag) {
|
||||
double ecoul = -prefactor * erfc_etar;
|
||||
force->pair->ev_tally(i, j, nlocal, newton_pair, 0., ecoul, 0., 0., 0., 0.);
|
||||
}
|
||||
if (vflag) v_tally(i, j, nlocal, newton_pair, fpair, delx, dely, delz);
|
||||
}
|
||||
}
|
||||
}
|
||||
@ -1621,3 +1625,70 @@ void FixElectrodeConp::unpack_forward_comm(int n, int first, double *buf)
|
||||
int const last = first + n;
|
||||
for (int i = first, m = 0; i < last; i++) atom->q[i] = buf[m++];
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Tally virial of pair interactions in pre_reverse. This cannot be done with pair->ev_tally()
|
||||
because compute_fdotr is called before pre_reverse, i.e. Virials need to be tallied even if fdotr
|
||||
is used.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void FixElectrodeConp::v_tally(int i, int j, int nlocal, int newton_pair, double fpair, double delx,
|
||||
double dely, double delz)
|
||||
{
|
||||
double v[6];
|
||||
if (vflag_either) {
|
||||
v[0] = delx * delx * fpair;
|
||||
v[1] = dely * dely * fpair;
|
||||
v[2] = delz * delz * fpair;
|
||||
v[3] = delx * dely * fpair;
|
||||
v[4] = delx * delz * fpair;
|
||||
v[5] = dely * delz * fpair;
|
||||
|
||||
if (vflag_global) {
|
||||
if (newton_pair) {
|
||||
virial[0] += v[0];
|
||||
virial[1] += v[1];
|
||||
virial[2] += v[2];
|
||||
virial[3] += v[3];
|
||||
virial[4] += v[4];
|
||||
virial[5] += v[5];
|
||||
} else {
|
||||
if (i < nlocal) {
|
||||
virial[0] += 0.5 * v[0];
|
||||
virial[1] += 0.5 * v[1];
|
||||
virial[2] += 0.5 * v[2];
|
||||
virial[3] += 0.5 * v[3];
|
||||
virial[4] += 0.5 * v[4];
|
||||
virial[5] += 0.5 * v[5];
|
||||
}
|
||||
if (j < nlocal) {
|
||||
virial[0] += 0.5 * v[0];
|
||||
virial[1] += 0.5 * v[1];
|
||||
virial[2] += 0.5 * v[2];
|
||||
virial[3] += 0.5 * v[3];
|
||||
virial[4] += 0.5 * v[4];
|
||||
virial[5] += 0.5 * v[5];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (vflag_atom) {
|
||||
if (newton_pair || i < nlocal) {
|
||||
vatom[i][0] += 0.5 * v[0];
|
||||
vatom[i][1] += 0.5 * v[1];
|
||||
vatom[i][2] += 0.5 * v[2];
|
||||
vatom[i][3] += 0.5 * v[3];
|
||||
vatom[i][4] += 0.5 * v[4];
|
||||
vatom[i][5] += 0.5 * v[5];
|
||||
}
|
||||
if (newton_pair || j < nlocal) {
|
||||
vatom[j][0] += 0.5 * v[0];
|
||||
vatom[j][1] += 0.5 * v[1];
|
||||
vatom[j][2] += 0.5 * v[2];
|
||||
vatom[j][3] += 0.5 * v[3];
|
||||
vatom[j][4] += 0.5 * v[4];
|
||||
vatom[j][5] += 0.5 * v[5];
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
@ -119,10 +119,11 @@ class FixElectrodeConp : public Fix {
|
||||
void create_taglist();
|
||||
void invert();
|
||||
void symmetrize();
|
||||
double gausscorr(int, bool);
|
||||
double gausscorr(int, int, bool);
|
||||
void update_charges();
|
||||
double potential_energy();
|
||||
double self_energy(int);
|
||||
void v_tally(int, int, int, int, double, double, double, double);
|
||||
void write_to_file(FILE *, const std::vector<tagint> &, const std::vector<std::vector<double>> &);
|
||||
void read_from_file(const std::string &input_file, double **, const std::string &);
|
||||
void compute_sd_vectors();
|
||||
|
||||
@ -32,6 +32,7 @@ ComputeMSDChunk::ComputeMSDChunk(LAMMPS *lmp, int narg, char **arg) :
|
||||
{
|
||||
if (narg != 4) error->all(FLERR, "Illegal compute msd/chunk command");
|
||||
|
||||
msdnchunk = 0;
|
||||
array_flag = 1;
|
||||
size_array_cols = 4;
|
||||
size_array_rows = 0;
|
||||
@ -117,14 +118,14 @@ void ComputeMSDChunk::compute_array()
|
||||
double massone;
|
||||
double unwrap[3];
|
||||
|
||||
int oldnchunk = nchunk;
|
||||
ComputeChunk::compute_array();
|
||||
int *ichunk = cchunk->ichunk;
|
||||
|
||||
// first time call, allocate per-chunk arrays
|
||||
// thereafter, require nchunk remain the same
|
||||
|
||||
if (!firstflag && (oldnchunk != nchunk))
|
||||
if (firstflag) msdnchunk = nchunk;
|
||||
else if (msdnchunk != nchunk)
|
||||
error->all(FLERR, "Compute msd/chunk nchunk is not static");
|
||||
|
||||
// zero local per-chunk values
|
||||
|
||||
@ -35,6 +35,8 @@ class ComputeMSDChunk : public ComputeChunk {
|
||||
double memory_usage() override;
|
||||
|
||||
private:
|
||||
int msdnchunk;
|
||||
|
||||
char *id_fix;
|
||||
class FixStoreGlobal *fix;
|
||||
|
||||
|
||||
@ -674,7 +674,9 @@ int DumpAtom::convert_image(int n, double *mybuf)
|
||||
memory->grow(sbuf,maxsbuf,"dump:sbuf");
|
||||
}
|
||||
|
||||
offset += sprintf(&sbuf[offset],format,
|
||||
offset += snprintf(&sbuf[offset],
|
||||
maxsbuf - offset,
|
||||
format,
|
||||
static_cast<tagint> (mybuf[m]),
|
||||
static_cast<int> (mybuf[m+1]),
|
||||
mybuf[m+2],mybuf[m+3],mybuf[m+4],
|
||||
@ -700,7 +702,9 @@ int DumpAtom::convert_noimage(int n, double *mybuf)
|
||||
memory->grow(sbuf,maxsbuf,"dump:sbuf");
|
||||
}
|
||||
|
||||
offset += sprintf(&sbuf[offset],format,
|
||||
offset += snprintf(&sbuf[offset],
|
||||
maxsbuf - offset,
|
||||
format,
|
||||
static_cast<tagint> (mybuf[m]),
|
||||
static_cast<int> (mybuf[m+1]),
|
||||
mybuf[m+2],mybuf[m+3],mybuf[m+4]);
|
||||
|
||||
@ -166,23 +166,24 @@ int DumpCFG::convert_string(int n, double *mybuf)
|
||||
}
|
||||
|
||||
for (j = 0; j < size_one; j++) {
|
||||
const auto maxsize = maxsbuf - offset;
|
||||
if (j == 0) {
|
||||
offset += sprintf(&sbuf[offset],"%f \n",mybuf[m]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,"%f \n",mybuf[m]);
|
||||
} else if (j == 1) {
|
||||
offset += sprintf(&sbuf[offset],"%s \n",typenames[(int) mybuf[m]]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,"%s \n",typenames[(int) mybuf[m]]);
|
||||
} else if (j >= 2) {
|
||||
if (vtype[j] == Dump::INT)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],static_cast<int> (mybuf[m]));
|
||||
else if (vtype[j] == Dump::DOUBLE)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],mybuf[m]);
|
||||
else if (vtype[j] == Dump::STRING)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],typenames[(int) mybuf[m]]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],typenames[(int) mybuf[m]]);
|
||||
else if (vtype[j] == Dump::BIGINT)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],static_cast<bigint> (mybuf[m]));
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],static_cast<bigint> (mybuf[m]));
|
||||
}
|
||||
m++;
|
||||
}
|
||||
offset += sprintf(&sbuf[offset],"\n");
|
||||
offset += snprintf(&sbuf[offset],maxsbuf-offset,"\n");
|
||||
}
|
||||
|
||||
} else if (unwrapflag == 1) {
|
||||
@ -195,29 +196,30 @@ int DumpCFG::convert_string(int n, double *mybuf)
|
||||
}
|
||||
|
||||
for (j = 0; j < size_one; j++) {
|
||||
const auto maxsize = maxsbuf - offset;
|
||||
if (j == 0) {
|
||||
offset += sprintf(&sbuf[offset],"%f \n",mybuf[m]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,"%f \n",mybuf[m]);
|
||||
} else if (j == 1) {
|
||||
offset += sprintf(&sbuf[offset],"%s \n",typenames[(int) mybuf[m]]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,"%s \n",typenames[(int) mybuf[m]]);
|
||||
} else if (j >= 2 && j <= 4) {
|
||||
unwrap_coord = (mybuf[m] - 0.5)/UNWRAPEXPAND + 0.5;
|
||||
offset += sprintf(&sbuf[offset],vformat[j],unwrap_coord);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],unwrap_coord);
|
||||
} else if (j >= 5) {
|
||||
if (vtype[j] == Dump::INT)
|
||||
offset +=
|
||||
sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
|
||||
snprintf(&sbuf[offset],maxsize,vformat[j],static_cast<int> (mybuf[m]));
|
||||
else if (vtype[j] == Dump::DOUBLE)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],mybuf[m]);
|
||||
else if (vtype[j] == Dump::STRING)
|
||||
offset +=
|
||||
sprintf(&sbuf[offset],vformat[j],typenames[(int) mybuf[m]]);
|
||||
snprintf(&sbuf[offset],maxsize,vformat[j],typenames[(int) mybuf[m]]);
|
||||
else if (vtype[j] == Dump::BIGINT)
|
||||
offset +=
|
||||
sprintf(&sbuf[offset],vformat[j],static_cast<bigint> (mybuf[m]));
|
||||
snprintf(&sbuf[offset],maxsize,vformat[j],static_cast<bigint> (mybuf[m]));
|
||||
}
|
||||
m++;
|
||||
}
|
||||
offset += sprintf(&sbuf[offset],"\n");
|
||||
offset += snprintf(&sbuf[offset],maxsbuf - offset,"\n");
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
@ -1358,20 +1358,21 @@ int DumpCustom::convert_string(int n, double *mybuf)
|
||||
}
|
||||
|
||||
for (j = 0; j < nfield; j++) {
|
||||
const auto maxsize = maxsbuf - offset;
|
||||
if (vtype[j] == Dump::INT)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],static_cast<int> (mybuf[m]));
|
||||
else if (vtype[j] == Dump::DOUBLE)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],mybuf[m]);
|
||||
else if (vtype[j] == Dump::STRING)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],typenames[(int) mybuf[m]]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],typenames[(int) mybuf[m]]);
|
||||
else if (vtype[j] == Dump::STRING2)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],atom->lmap->typelabel[(int) mybuf[m]-1].c_str());
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],atom->lmap->typelabel[(int) mybuf[m]-1].c_str());
|
||||
else if (vtype[j] == Dump::BIGINT)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],
|
||||
static_cast<bigint> (mybuf[m]));
|
||||
m++;
|
||||
}
|
||||
offset += sprintf(&sbuf[offset],"\n");
|
||||
offset += snprintf(&sbuf[offset],maxsbuf-offset,"\n");
|
||||
}
|
||||
|
||||
return offset;
|
||||
@ -1909,9 +1910,9 @@ int DumpCustom::modify_param(int narg, char **arg)
|
||||
if (ptr == nullptr)
|
||||
error->all(FLERR,"Dump_modify int format does not contain d character");
|
||||
char str[8];
|
||||
sprintf(str,"%s",BIGINT_FORMAT);
|
||||
snprintf(str,8,"%s",BIGINT_FORMAT);
|
||||
*ptr = '\0';
|
||||
sprintf(format_bigint_user,"%s%s%s",format_int_user,&str[1],ptr+1);
|
||||
snprintf(format_bigint_user,n,"%s%s%s",format_int_user,&str[1],ptr+1);
|
||||
*ptr = 'd';
|
||||
|
||||
} else if (strcmp(arg[1],"float") == 0) {
|
||||
|
||||
@ -590,15 +590,16 @@ int DumpGrid::convert_string(int n, double *mybuf)
|
||||
}
|
||||
|
||||
for (j = 0; j < nfield; j++) {
|
||||
const auto maxsize = maxsbuf - offset;
|
||||
if (vtype[j] == Dump::INT)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],static_cast<int> (mybuf[m]));
|
||||
else if (vtype[j] == Dump::DOUBLE)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],mybuf[m]);
|
||||
else if (vtype[j] == Dump::BIGINT)
|
||||
offset += sprintf(&sbuf[offset],vformat[j], static_cast<bigint> (mybuf[m]));
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j], static_cast<bigint> (mybuf[m]));
|
||||
m++;
|
||||
}
|
||||
offset += sprintf(&sbuf[offset],"\n");
|
||||
offset += snprintf(&sbuf[offset],maxsbuf-offset,"\n");
|
||||
}
|
||||
|
||||
return offset;
|
||||
@ -776,9 +777,9 @@ int DumpGrid::modify_param(int narg, char **arg)
|
||||
if (ptr == nullptr)
|
||||
error->all(FLERR,"Dump_modify int format does not contain d character");
|
||||
char str[8];
|
||||
sprintf(str,"%s",BIGINT_FORMAT);
|
||||
snprintf(str,8,"%s",BIGINT_FORMAT);
|
||||
*ptr = '\0';
|
||||
sprintf(format_bigint_user,"%s%s%s",format_int_user,&str[1],ptr+1);
|
||||
snprintf(format_bigint_user,n,"%s%s%s",format_int_user,&str[1],ptr+1);
|
||||
*ptr = 'd';
|
||||
|
||||
} else if (strcmp(arg[1],"float") == 0) {
|
||||
|
||||
@ -264,9 +264,9 @@ int DumpLocal::modify_param(int narg, char **arg)
|
||||
if (ptr == nullptr)
|
||||
error->all(FLERR, "Dump_modify int format does not contain d character");
|
||||
char str[8];
|
||||
sprintf(str,"%s",BIGINT_FORMAT);
|
||||
snprintf(str,8,"%s",BIGINT_FORMAT);
|
||||
*ptr = '\0';
|
||||
sprintf(format_bigint_user,"%s%s%s",format_int_user,&str[1],ptr+1);
|
||||
snprintf(format_bigint_user,n,"%s%s%s",format_int_user,&str[1],ptr+1);
|
||||
*ptr = 'd';
|
||||
|
||||
} else if (strcmp(arg[1],"float") == 0) {
|
||||
@ -387,17 +387,18 @@ int DumpLocal::convert_string(int n, double *mybuf)
|
||||
}
|
||||
|
||||
for (j = 0; j < size_one; j++) {
|
||||
const auto maxsize = maxsbuf - offset;
|
||||
if (vtype[j] == Dump::INT)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],static_cast<int> (mybuf[m]));
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],static_cast<int> (mybuf[m]));
|
||||
else if (vtype[j] == Dump::DOUBLE)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],mybuf[m]);
|
||||
else if (vtype[j] == Dump::BIGINT)
|
||||
offset += sprintf(&sbuf[offset],vformat[j],static_cast<bigint> (mybuf[m]));
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],static_cast<bigint> (mybuf[m]));
|
||||
else
|
||||
offset += sprintf(&sbuf[offset],vformat[j],mybuf[m]);
|
||||
offset += snprintf(&sbuf[offset],maxsize,vformat[j],mybuf[m]);
|
||||
m++;
|
||||
}
|
||||
offset += sprintf(&sbuf[offset],"\n");
|
||||
offset += snprintf(&sbuf[offset],maxsbuf-offset,"\n");
|
||||
}
|
||||
|
||||
return offset;
|
||||
|
||||
@ -85,7 +85,7 @@ void DumpXYZ::init_style()
|
||||
typenames = new char*[ntypes+1];
|
||||
for (int itype = 1; itype <= ntypes; itype++) {
|
||||
typenames[itype] = new char[12];
|
||||
sprintf(typenames[itype],"%d",itype);
|
||||
snprintf(typenames[itype],12,"%d",itype);
|
||||
}
|
||||
}
|
||||
|
||||
@ -206,7 +206,7 @@ int DumpXYZ::convert_string(int n, double *mybuf)
|
||||
memory->grow(sbuf,maxsbuf,"dump:sbuf");
|
||||
}
|
||||
|
||||
offset += sprintf(&sbuf[offset], format, typenames[static_cast<int> (mybuf[m+1])],
|
||||
offset += snprintf(&sbuf[offset], maxsbuf-offset, format, typenames[static_cast<int> (mybuf[m+1])],
|
||||
mybuf[m+2], mybuf[m+3], mybuf[m+4]);
|
||||
m += size_one;
|
||||
}
|
||||
|
||||
@ -1025,7 +1025,7 @@ char *Variable::retrieve(const char *name)
|
||||
error->all(FLERR, "Variable {}: format variable {} has incompatible style",
|
||||
names[ivar],data[ivar][0]);
|
||||
double answer = compute_equal(jvar);
|
||||
sprintf(data[ivar][2],data[ivar][1],answer);
|
||||
snprintf(data[ivar][2],VALUELENGTH,data[ivar][1],answer);
|
||||
str = data[ivar][2];
|
||||
|
||||
} else if (style[ivar] == GETENV) {
|
||||
|
||||
@ -121,7 +121,7 @@ TEST_F(LAMMPS_extract_variable, loop_pad)
|
||||
char str[10];
|
||||
char *fstr;
|
||||
for (i = 1; i <= 10; i++) {
|
||||
std::sprintf(str, "%02d", i);
|
||||
std::snprintf(str, 10, "%02d", i);
|
||||
fstr = f_lammps_extract_variable_loop_pad();
|
||||
EXPECT_STREQ(fstr, str);
|
||||
std::free(fstr);
|
||||
@ -170,7 +170,7 @@ TEST_F(LAMMPS_extract_variable, format)
|
||||
char str[16];
|
||||
char *fstr;
|
||||
for (i = 1; i <= 10; i++) {
|
||||
std::sprintf(str, "%.6G", std::exp(i));
|
||||
std::snprintf(str, 16, "%.6G", std::exp(i));
|
||||
fstr = f_lammps_extract_variable_format();
|
||||
EXPECT_STREQ(fstr, str);
|
||||
std::free(fstr);
|
||||
@ -185,7 +185,7 @@ TEST_F(LAMMPS_extract_variable, format_pad)
|
||||
char str[16];
|
||||
char *fstr;
|
||||
for (i = 1; i <= 10; i++) {
|
||||
std::sprintf(str, "%08.6G", std::exp(i));
|
||||
std::snprintf(str, 16, "%08.6G", std::exp(i));
|
||||
fstr = f_lammps_extract_variable_format_pad();
|
||||
EXPECT_STREQ(fstr, str);
|
||||
std::free(fstr);
|
||||
|
||||
Reference in New Issue
Block a user