Merge pull request #671 from danicholson/user-uef
USER-UEF package for extensional flow
This commit is contained in:
Steve Plimpton
@ -105,7 +105,7 @@ set(OTHER_PACKAGES KIM PYTHON MSCG MPIIO VORONOI POEMS LATTE
|
||||
USER-CGSDK USER-COLVARS USER-DIFFRACTION USER-DPD USER-DRUDE USER-EFF
|
||||
USER-FEP USER-H5MD USER-LB USER-MANIFOLD USER-MEAMC USER-MGPT USER-MISC
|
||||
USER-MOLFILE USER-NETCDF USER-PHONON USER-QTB USER-REAXC USER-SMD
|
||||
USER-SMTBQ USER-SPH USER-TALLY USER-VTK USER-QUIP USER-QMMM)
|
||||
USER-SMTBQ USER-SPH USER-TALLY USER-UEF USER-VTK USER-QUIP USER-QMMM)
|
||||
set(ACCEL_PACKAGES USER-OMP KOKKOS OPT USER-INTEL GPU)
|
||||
foreach(PKG ${DEFAULT_PACKAGES})
|
||||
option(ENABLE_${PKG} "Build ${PKG} Package" ${ENABLE_ALL})
|
||||
|
||||
BIN
doc/src/JPG/uef_frames.jpg
Normal file
BIN
doc/src/JPG/uef_frames.jpg
Normal file
Binary file not shown.
|
After Width: | Height: | Size: 12 KiB |
@ -720,6 +720,8 @@ package"_Section_start.html#start_3.
|
||||
"nve/eff"_fix_nve_eff.html,
|
||||
"nvt/eff"_fix_nh_eff.html,
|
||||
"nvt/sllod/eff"_fix_nvt_sllod_eff.html,
|
||||
"npt/uef"_fix_nh_uef.html,
|
||||
"nvt/uef"_fix_nh_uef.html,
|
||||
"phonon"_fix_phonon.html,
|
||||
"pimd"_fix_pimd.html,
|
||||
"qbmsst"_fix_qbmsst.html,
|
||||
@ -857,6 +859,7 @@ package"_Section_start.html#start_3.
|
||||
"meso/t/atom"_compute_meso_t_atom.html,
|
||||
"pe/tally"_compute_tally.html,
|
||||
"pe/mol/tally"_compute_tally.html,
|
||||
"pressure/uef"_compute_pressure_uef.html,
|
||||
"saed"_compute_saed.html,
|
||||
"smd/contact/radius"_compute_smd_contact_radius.html,
|
||||
"smd/damage"_compute_smd_damage.html,
|
||||
@ -885,6 +888,7 @@ package"_Section_start.html#start_3.
|
||||
"temp/deform/eff"_compute_temp_deform_eff.html,
|
||||
"temp/region/eff"_compute_temp_region_eff.html,
|
||||
"temp/rotate"_compute_temp_rotate.html,
|
||||
"temp/uef"_compute_temp_uef.html,
|
||||
"xrd"_compute_xrd.html :tb(c=6,ea=c)
|
||||
|
||||
:line
|
||||
|
||||
@ -150,6 +150,7 @@ Package, Description, Doc page, Example, Library
|
||||
"USER-SMTBQ"_#USER-SMTBQ, second moment tight binding QEq potential,"pair_style smtbq"_pair_smtbq.html, USER/smtbq, -
|
||||
"USER-SPH"_#USER-SPH, smoothed particle hydrodynamics,"SPH User Guide"_PDF/SPH_LAMMPS_userguide.pdf, USER/sph, -
|
||||
"USER-TALLY"_#USER-TALLY, pairwise tally computes,"compute XXX/tally"_compute_tally.html, USER/tally, -
|
||||
"USER-UEF"_#USER-UEF, extensional flow,"fix nvt/uef"_fix_nh_uef.html, USER/uef, -
|
||||
"USER-VTK"_#USER-VTK, dump output via VTK, "compute vtk"_dump_vtk.html, -, ext :tb(ea=c,ca1=l)
|
||||
|
||||
:line
|
||||
@ -2771,6 +2772,37 @@ examples/USER/tally :ul
|
||||
|
||||
:line
|
||||
|
||||
USER-UEF package :link(USER-UEF),h4
|
||||
|
||||
[Contents:]
|
||||
|
||||
A fix style for the integration of the equations of motion under
|
||||
extensional flow with proper boundary conditions, as well
|
||||
as several supporting compute styles and an output option.
|
||||
|
||||
[Author:] David Nicholson (MIT).
|
||||
|
||||
[Install or un-install:]
|
||||
|
||||
make yes-user-uef
|
||||
make machine :pre
|
||||
|
||||
make no-user-uef
|
||||
make machine :pre
|
||||
|
||||
[Supporting info:]
|
||||
|
||||
src/USER-UEF: filenames -> commands
|
||||
src/USER-UEF/README
|
||||
"fix nvt/uef"_fix_nh_uef.html
|
||||
"fix npt/uef"_fix_nh_uef.html
|
||||
"compute pressure/uef"_compute_pressure_uef.html
|
||||
"compute temp/uef"_compute_temp_uef.html
|
||||
"dump cfg/uef"_dump_cfg_uef.html
|
||||
examples/uef :ul
|
||||
|
||||
:line
|
||||
|
||||
USER-VTK package :link(USER-VTK),h4
|
||||
|
||||
[Contents:]
|
||||
|
||||
@ -32,6 +32,7 @@ Commands :h1
|
||||
dimension
|
||||
displace_atoms
|
||||
dump
|
||||
dump_cfg_uef
|
||||
dump_h5md
|
||||
dump_image
|
||||
dump_modify
|
||||
|
||||
61
doc/src/compute_pressure_uef.txt
Normal file
61
doc/src/compute_pressure_uef.txt
Normal file
@ -0,0 +1,61 @@
|
||||
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
|
||||
|
||||
:link(lws,http://lammps.sandia.gov)
|
||||
:link(ld,Manual.html)
|
||||
:link(lc,Section_commands.html#comm)
|
||||
|
||||
:line
|
||||
|
||||
compute pressure/uef command :h3
|
||||
|
||||
[Syntax:]
|
||||
|
||||
compute ID group-ID pressure/uef temp-ID keyword ... :pre
|
||||
|
||||
ID, group-ID are documented in "compute"_compute.html command
|
||||
pressure/uef = style name of this compute command
|
||||
temp-ID = ID of compute that calculates temperature, can be NULL if not needed
|
||||
zero or more keywords may be appended
|
||||
keyword = {ke} or {pair} or {bond} or {angle} or {dihedral} or {improper} or {kspace} or {fix} or {virial} :ul
|
||||
|
||||
[Examples:]
|
||||
|
||||
compute 1 all pressure/uef my_temp_uef
|
||||
compute 2 all pressure/uef my_temp_uef virial :pre
|
||||
|
||||
[Description:]
|
||||
|
||||
This command is used to compute the pressure tensor in
|
||||
the reference frame of the applied flow field when
|
||||
"fix nvt/uef"_fix_nh_uef.html" or
|
||||
"fix npt/uef"_fix_nh_uef.html" is used.
|
||||
It is not necessary to use this command to compute the scalar
|
||||
value of the pressure. A "compute pressure"_compute_pressure.html
|
||||
may be used for that purpose.
|
||||
|
||||
The keywords and output information are documented in
|
||||
"compute_pressure"_compute_pressure.html.
|
||||
|
||||
[Restrictions:]
|
||||
|
||||
This fix is part of the USER-UEF package. It is only enabled if
|
||||
LAMMPS was built with that package. See the
|
||||
"Making LAMMPS"_Section_start.html#start_3 section for more info.
|
||||
|
||||
This command can only be used when "fix nvt/uef"_fix_nh_uef.html
|
||||
or "fix npt/uef"_fix_nh_uef.html is active.
|
||||
|
||||
The kinetic contribution to the pressure tensor
|
||||
will be accurate only when
|
||||
the compute specificed by {temp-ID} is a
|
||||
"compute temp/uef"_compute_temp_uef.html.
|
||||
|
||||
[Related commands:]
|
||||
|
||||
"compute pressure"_compute_pressure.html,
|
||||
"fix nvt/uef"_fix_nh_uef.html,
|
||||
"compute temp/uef"_compute_temp_uef.html
|
||||
|
||||
[Default:] none
|
||||
|
||||
|
||||
52
doc/src/compute_temp_uef.txt
Normal file
52
doc/src/compute_temp_uef.txt
Normal file
@ -0,0 +1,52 @@
|
||||
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
|
||||
|
||||
:link(lws,http://lammps.sandia.gov)
|
||||
:link(ld,Manual.html)
|
||||
:link(lc,Section_commands.html#comm)
|
||||
|
||||
:line
|
||||
|
||||
compute temp/uef command :h3
|
||||
|
||||
[Syntax:]
|
||||
|
||||
compute ID group-ID temp/uef :pre
|
||||
|
||||
ID, group-ID are documented in "compute"_compute.html command
|
||||
temp/uef = style name of this compute command :ul
|
||||
|
||||
[Examples:]
|
||||
|
||||
compute 1 all temp/uef
|
||||
compute 2 sel temp/uef :pre
|
||||
|
||||
[Description:]
|
||||
|
||||
This command is used to compute the kinetic energy tensor in
|
||||
the reference frame of the applied flow field when
|
||||
"fix nvt/uef"_fix_nh_uef.html" or
|
||||
"fix npt/uef"_fix_nh_uef.html" is used.
|
||||
It is not necessary to use this command to compute the scalar
|
||||
value of the temperature. A "compute temp"_compute_temp.html
|
||||
may be used for that purpose.
|
||||
|
||||
Output information for this command can be found in the
|
||||
documentation for "compute temp"_compute_temp.html.
|
||||
|
||||
[Restrictions:]
|
||||
|
||||
This fix is part of the USER-UEF package. It is only enabled if
|
||||
LAMMPS was built with that package. See the
|
||||
"Making LAMMPS"_Section_start.html#start_3 section for more info.
|
||||
|
||||
This command can only be used when "fix nvt/uef"_fix_nh_uef.html
|
||||
or "fix npt/uef"_fix_nh_uef.html is active.
|
||||
|
||||
[Related commands:]
|
||||
|
||||
"compute temp"_compute_temp.html,
|
||||
"fix nvt/uef"_fix_nh_uef.html,
|
||||
"compute pressure/uef"_compute_pressure_uef.html
|
||||
|
||||
|
||||
[Default:] none
|
||||
@ -65,6 +65,7 @@ Computes :h1
|
||||
compute_pe_atom
|
||||
compute_plasticity_atom
|
||||
compute_pressure
|
||||
compute_pressure_uef
|
||||
compute_property_atom
|
||||
compute_property_chunk
|
||||
compute_property_local
|
||||
@ -114,6 +115,7 @@ Computes :h1
|
||||
compute_temp_region_eff
|
||||
compute_temp_rotate
|
||||
compute_temp_sphere
|
||||
compute_temp_uef
|
||||
compute_ti
|
||||
compute_torque_chunk
|
||||
compute_vacf
|
||||
|
||||
53
doc/src/dump_cfg_uef.txt
Normal file
53
doc/src/dump_cfg_uef.txt
Normal file
@ -0,0 +1,53 @@
|
||||
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
|
||||
|
||||
:link(lws,http://lammps.sandia.gov)
|
||||
:link(ld,Manual.html)
|
||||
:link(lc,Section_commands.html#comm)
|
||||
|
||||
:line
|
||||
|
||||
dump cfg/uef command :h3
|
||||
|
||||
[Syntax:]
|
||||
|
||||
dump ID group-ID cfg/uef N file mass type xs ys zs args :pre
|
||||
|
||||
ID = user-assigned name for the dump :ulb,l
|
||||
group-ID = ID of the group of atoms to be dumped :l
|
||||
N = dump every this many timesteps :l
|
||||
file = name of file to write dump info to :l
|
||||
args = same as args for "dump custom"_dump.html :pre
|
||||
|
||||
:ule
|
||||
|
||||
[Examples:]
|
||||
|
||||
dump 1 all cfg/uef 10 dump.*.cfg mass type xs ys zs
|
||||
dump 2 all cfg/uef 100 dump.*.cfg mass type xs ys zs id c_stress :pre
|
||||
|
||||
[Description:]
|
||||
|
||||
This command is used to dump atomic coordinates in the
|
||||
reference frame of the applied flow field when
|
||||
"fix nvt/uef"_fix_nh_uef.html or
|
||||
"fix npt/uef"_fix_nh_uef.html or is used. Only the atomic
|
||||
coordinates and frame-invariant scalar quantities
|
||||
will be in the flow frame. If velocities are selected
|
||||
as output, for example, they will not be in the same
|
||||
reference frame as the atomic positions.
|
||||
|
||||
[Restrictions:]
|
||||
|
||||
This fix is part of the USER-UEF package. It is only enabled if
|
||||
LAMMPS was built with that package. See the
|
||||
"Making LAMMPS"_Section_start.html#start_3 section for more info.
|
||||
|
||||
This command can only be used when "fix nvt/uef"_fix_nh_uef.html
|
||||
or "fix npt/uef"_fix_nh_uef.html is active.
|
||||
|
||||
[Related commands:]
|
||||
|
||||
"dump"_dump.html,
|
||||
"fix nvt/uef"_fix_nh_uef.html
|
||||
|
||||
[Default:] none
|
||||
211
doc/src/fix_nh_uef.txt
Normal file
211
doc/src/fix_nh_uef.txt
Normal file
@ -0,0 +1,211 @@
|
||||
<"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
|
||||
|
||||
:link(lws,http://lammps.sandia.gov)
|
||||
:link(ld,Manual.html)
|
||||
:link(lc,Section_commands.html#comm)
|
||||
|
||||
:line
|
||||
|
||||
fix nvt/uef command :h3
|
||||
fix npt/uef command :h3
|
||||
|
||||
[Syntax:]
|
||||
|
||||
fix ID group-ID style_name erate edot_x edot_y temp Tstart Tstop Tdamp keyword value ... :pre
|
||||
|
||||
ID, group-ID are documented in "fix"_fix.html command :ulb,l
|
||||
style_name = {nvt/uef} or {npt/uef} :l
|
||||
{Tstart}, {Tstop}, and {Tdamp} are documented in the "fix npt"_fix_nh.html command :l
|
||||
{edot_x} and {edot_y} are the strain rates in the x and y directions (1/(time units)) :l
|
||||
one or more keyword/value pairs may be appended :l
|
||||
keyword = {ext} or {strain} or {iso} or {x} or {y} or {z} or {tchain} or {pchain} or {tloop} or {ploop} or {mtk}
|
||||
{ext} value = {x} or {y} or {z} or {xy} or {yz} or {xz} = external dimensions
|
||||
This keyword sets the external dimensions used to calculate the scalar pressure
|
||||
{strain} values = e_x e_y = initial strain
|
||||
Use of this keyword is usually not necessary, but may be needed to resume a run with a data file.
|
||||
The {iso}, {x}, {y}, {z}, {tchain}, {pchain}, {tloop}, {ploop}, and {mtk} keywords are documented in the "fix npt"_fix_nh.html command. :pre
|
||||
:ule
|
||||
|
||||
[Examples:]
|
||||
|
||||
fix uniax_nvt all nvt/uef temp 400 400 100 erate 0.00001 -0.000005
|
||||
fix biax_nvt all nvt/uef temp 400 400 100 erate 0.000005 0.000005
|
||||
fix uniax_npt all npt/uef temp 400 400 300 iso 1 1 3000 erate 0.00001 -0.000005 ext yz
|
||||
fix biax_npt all npt/uef temp 400 400 100 erate -0.00001 0.000005 x 1 1 3000 :pre
|
||||
|
||||
[Description:]
|
||||
|
||||
|
||||
This fix is used to simulate non-equilibrium molecular dynamics (NEMD)
|
||||
under diagonal flow fields, including uniaxial and biaxial flow.
|
||||
Simulations under extensional flow may be carried out for an
|
||||
indefinite amount of time. It is an implementation of the boundary
|
||||
conditions from "(Dobson)"_#Dobson, and also uses numerical lattice reduction
|
||||
as was proposed by "(Hunt)"_#Hunt. The lattice reduction algorithm is from
|
||||
"(Semaev)"_Semaev. The fix is intended for simulations of
|
||||
homogeneous flows, and integrates the SLLOD equations of motion,
|
||||
originally proposed by Hoover and Ladd (see "(Evans and Morriss)"_#Sllod).
|
||||
Additional detail about this implementation can be found in
|
||||
"(Nicholson and Rutledge)"_#Nicholson.
|
||||
|
||||
The applied flow field is set by the {eps} keyword. The values {edot_x}
|
||||
and {edot_y} correspond to the strain rates in the xx and yy directions.
|
||||
It is implicitly assumed that the flow field is traceless, and therefore
|
||||
the strain rate in the zz direction is eqal to -({edot_x} + {edot_y}).
|
||||
|
||||
NOTE: Due to an instability in the SLLOD equations under extension,
|
||||
"fix momentum"_fix_momentum.html should be used to regularly reset the
|
||||
linear momentum.
|
||||
|
||||
The boundary conditions require a simulation box that does not have a
|
||||
consistent alignment relative to the applied flow field. Since LAMMPS
|
||||
utilizes an upper-triangular simulation box, it is not possible to express
|
||||
the evolving simulation box in the same coordinate system as the flow field.
|
||||
This fix keeps track of two coordinate systems: the flow frame, and the
|
||||
upper triangular LAMMPS frame. The coordinate systems are related to
|
||||
each other through the QR decomposition, as is illustrated in the image below.
|
||||
|
||||
:c,image(JPG/uef_frames.jpg)
|
||||
|
||||
During most molecular dynamics operations, the system is represented in the
|
||||
LAMMPS frame. Only when the positions and velocities are updated is the system
|
||||
rotated to the flow frame, and it is rotated back to the LAMMPS frame
|
||||
immediately afterwards. For this reason, all vector-valued quantities
|
||||
(except for the tensors from "compute_pressure/uef"_compute_pressure_uef.html
|
||||
and "compute_temp/uef"_compute_temp_uef.html) will
|
||||
be computed in the LAMMPS frame. Rotationally invariant scalar quantities like
|
||||
the temperature and hydrostatic pressure are frame-invariant and will be
|
||||
computed correctly. Additionally, the system is in the LAMMPS frame during all of the
|
||||
output steps, and therefore trajectory files made using the dump command
|
||||
will be in the LAMMPS frame unless the "dump_cfg/uef"_dump_cfg_uef.html command is used.
|
||||
|
||||
:line
|
||||
|
||||
Temperature control is achieved with the default Nose-Hoover style
|
||||
thermostat documented in "fix npt"_fix_nh.html. When this fix is active,
|
||||
only the peculiar velocity of each atom is stored, defined as the velocity
|
||||
relative to the streaming velocity. This is in contrast to
|
||||
"fix nvt/sllod"_fix_nvt_sllod.html, which uses a lab-frame velocity, and
|
||||
removes the contribution from the streaming velocity in order to compute
|
||||
the temperature.
|
||||
|
||||
Pressure control is achieved using the default Nose-Hoover barostat documented
|
||||
in "fix npt"_fix_nh.html. There are two ways to control the pressure using this
|
||||
fix. The first method involves using the {ext} keyword along with the {iso} pressure
|
||||
style. With this method, the pressure is controlled by scaling the simulation box
|
||||
isotropically to achieve the average pressure only in the directions specified by {ext}.
|
||||
For example, if the {ext} value is set to {xy}, the average pressure (Pxx+Pyy)/2
|
||||
will be controlled.
|
||||
|
||||
This example command will control the total hydrostatic pressure under uniaxial tension:
|
||||
|
||||
fix f1 all npt/uef temp 0.7 0.7 0.5 iso 1 1 5 erate -0.5 -0.5 ext xyz :pre
|
||||
|
||||
This example command will control the average stress in compression directions, which would
|
||||
typically correspond to free surfaces under drawing with uniaxial tension:
|
||||
|
||||
fix f2 all npt/uef temp 0.7 0.7 0.5 iso 1 1 5 erate -0.5 -0.5 ext xy :pre
|
||||
|
||||
The second method for pressure control involves setting the normal stresses using
|
||||
the {x}, {y} , and/or {z} keywords. When using this method, the same pressure must be
|
||||
specified via {Pstart} and {Pstop} for all dimensions controlled. Any choice of
|
||||
pressure conditions that would cause LAMMPS to compute a deviatoric stress are not
|
||||
permissible and will result in an error. Additionally, all dimensions with
|
||||
controlled stress must have the same applied strain rate. The {ext} keyword must be
|
||||
set to the default value ({xyz}) when using this method.
|
||||
|
||||
For example, the following commands will work:
|
||||
|
||||
fix f3 all npt/uef temp 0.7 0.7 0.5 x 1 1 5 y 1 1 5 erate -0.5 -0.5
|
||||
fix f4 all npt/uef temp 0.7 0.7 0.5 z 1 1 5 erate 0.5 0.5 :pre
|
||||
|
||||
The following commands will not work:
|
||||
|
||||
fix f5 all npt/uef temp 0.7 0.7 0.5 x 1 1 5 z 1 1 5 erate -0.5 -0.5
|
||||
fix f6 all npt/uef temp 0.7 0.7 0.5 x 1 1 5 z 2 2 5 erate 0.5 0.5 :pre
|
||||
|
||||
:line
|
||||
|
||||
These fix computes a temperature and pressure each timestep. To do
|
||||
this, it creates its own computes of style "temp/uef" and "pressure/uef",
|
||||
as if one of these two sets of commands had been issued:
|
||||
|
||||
compute fix-ID_temp group-ID temp/uef
|
||||
compute fix-ID_press group-ID pressure/uef fix-ID_temp :pre
|
||||
|
||||
compute fix-ID_temp all temp/uef
|
||||
compute fix-ID_press all pressure/uef fix-ID_temp :pre
|
||||
|
||||
See the "compute temp/uef"_compute_temp_uef.html and "compute
|
||||
pressure/uef"_compute_pressure_uef.html commands for details. Note that the
|
||||
IDs of the new computes are the fix-ID + underscore + "temp" or fix_ID
|
||||
+ underscore + "press".
|
||||
|
||||
[Restart, fix_modify, output, run start/stop, minimize info:]
|
||||
|
||||
The fix writes the state of all the thermostat and barostat
|
||||
variables, as well as the cumulative strain applied, to
|
||||
"binary restart files"_restart.html. See the
|
||||
"read_restart"_read_restart.html command for info on how to re-specify
|
||||
a fix in an input script that reads a restart file, so that the
|
||||
operation of the fix continues in an uninterrupted fashion.
|
||||
|
||||
NOTE: It is not necessary to set the {strain} keyword when resuming
|
||||
a run from a restart file. Only for resuming from data files,
|
||||
which do not contain the cumulative applied strain, will
|
||||
this keyword be necessary.
|
||||
|
||||
This fix can be used with the "fix_modify"_fix_modify.html
|
||||
{temp} and {press} options. The temperature and pressure computes
|
||||
used must be of type {temp/uef} and {pressure/uef}.
|
||||
|
||||
This fix computes the same global scalar and vecor quantities
|
||||
as "fix npt"_fix_nh.html.
|
||||
|
||||
The fix is not invoked during "energy
|
||||
minimization"_minimize.html.
|
||||
|
||||
[Restrictions:]
|
||||
|
||||
This fix is part of the USER-UEF package. It is only enabled if
|
||||
LAMMPS was built with that package. See the
|
||||
"Making LAMMPS"_Section_start.html#start_3 section for more info.
|
||||
|
||||
Due to requirements of the boundary conditions, when the {strain} keyword
|
||||
is set to zero (or unset), the initial simulation box must be cubic and
|
||||
have style triclinic. If the box is initially of type ortho, use
|
||||
"change_box"_change_box.html before invoking the fix.
|
||||
|
||||
NOTE: When resuming from restart files, you may need to use "box tilt large"_box.html
|
||||
since lammps has internal criteria from lattice reduction that are not
|
||||
the same as the criteria in the numerical lattice reduction algorithm.
|
||||
|
||||
[Related commands:]
|
||||
"fix nvt"_fix_nh.html,
|
||||
"fix nvt/sllod"_fix_nvt_sllod.html,
|
||||
"compute temp/uef"_compute_temp_uef.html,
|
||||
"compute pressure/uef"_compute_pressure_uef.html,
|
||||
"dump cfg/uef"_dump_cfg_uef.html
|
||||
|
||||
[Default:]
|
||||
|
||||
The default keyword values specific to this fix are exy = xyz, strain = 0 0.
|
||||
The remaining defaults are the same as for {fix npt}_fix_nh.html except tchain = 1.
|
||||
The reason for this change is given in "fix nvt/sllod"_fix_nvt_sllod.html.
|
||||
|
||||
:line
|
||||
|
||||
:link(Dobson)
|
||||
[(Dobson)] Dobson, J Chem Phys, 141, 184103 (2014).
|
||||
|
||||
:link(Hunt)
|
||||
[(Hunt)] Hunt, Mol Simul, 42, 347 (2016).
|
||||
|
||||
:link(Semaev)
|
||||
[(Semaev)] Semaev, Cryptography and Lattices, 181 (2001).
|
||||
|
||||
:link(Sllod)
|
||||
[(Evans and Morriss)] Evans and Morriss, Phys Rev A, 30, 1528 (1984).
|
||||
|
||||
:link(Nicholson)
|
||||
[(Nicholson and Rutledge)] Nicholson and Rutledge, J Chem Phys, 145, 244903 (2016).
|
||||
@ -77,6 +77,7 @@ Fixes :h1
|
||||
fix_neb
|
||||
fix_nh
|
||||
fix_nh_eff
|
||||
fix_nh_uef
|
||||
fix_nph_asphere
|
||||
fix_nph_body
|
||||
fix_nph_sphere
|
||||
|
||||
@ -62,6 +62,7 @@ dump_modify.html
|
||||
dump_molfile.html
|
||||
dump_netcdf.html
|
||||
dump_vtk.html
|
||||
dump_cfg_uef.html
|
||||
echo.html
|
||||
fix.html
|
||||
fix_modify.html
|
||||
@ -232,6 +233,7 @@ fix_nvt_manifold_rattle.html
|
||||
fix_nvt_sllod.html
|
||||
fix_nvt_sllod_eff.html
|
||||
fix_nvt_sphere.html
|
||||
fix_nh_uef.html
|
||||
fix_oneway.html
|
||||
fix_orient.html
|
||||
fix_phonon.html
|
||||
@ -356,6 +358,7 @@ compute_pe.html
|
||||
compute_pe_atom.html
|
||||
compute_plasticity_atom.html
|
||||
compute_pressure.html
|
||||
compute_pressure_uef.html
|
||||
compute_property_atom.html
|
||||
compute_property_chunk.html
|
||||
compute_property_local.html
|
||||
@ -405,6 +408,7 @@ compute_temp_region.html
|
||||
compute_temp_region_eff.html
|
||||
compute_temp_rotate.html
|
||||
compute_temp_sphere.html
|
||||
compute_temp_uef.html
|
||||
compute_ti.html
|
||||
compute_torque_chunk.html
|
||||
compute_vacf.html
|
||||
|
||||
19
examples/USER/uef/README
Normal file
19
examples/USER/uef/README
Normal file
@ -0,0 +1,19 @@
|
||||
This directory contains two short example scripts for the USER-UEF
|
||||
package.
|
||||
|
||||
|
||||
nvt_uniaxial:
|
||||
A simulation of a WCA fluid under constant volume uniaxial strain
|
||||
with stretching in the z direction, compression in the x and y
|
||||
directions.The xx, yy and zz components of the pressure tensor are
|
||||
written to the log file. The simulation conditions are the same
|
||||
used in the paper by Thomas Hunt [http://arxiv.org/abs/1310.3905]
|
||||
and should yield an zz pressure near 5 and xx/yy pressures near 10.
|
||||
|
||||
npt_biaxial:
|
||||
A simulation of a WCA fluid under stress-controlled biaxial strain
|
||||
with stretching in the x and y directions, compression in the z
|
||||
direction. The zz stress is controlled using the ext keyword in fix
|
||||
npt/uef. The xx yy and zz components of the pressure tensor are
|
||||
written to the log file. The simulations conditions should yield
|
||||
xx/yy pressures near 5 and zz pressure near 10
|
||||
1022
examples/USER/uef/npt_biaxial/data.wca
Normal file
1022
examples/USER/uef/npt_biaxial/data.wca
Normal file
File diff suppressed because it is too large
Load Diff
31
examples/USER/uef/npt_biaxial/in.npt_biaxial
Normal file
31
examples/USER/uef/npt_biaxial/in.npt_biaxial
Normal file
@ -0,0 +1,31 @@
|
||||
# biaxial NPT deformation of WCA fluid
|
||||
|
||||
units lj
|
||||
atom_style atomic
|
||||
|
||||
|
||||
pair_style lj/cut 1.122562
|
||||
read_data data.wca
|
||||
pair_coeff 1 1 1.0 1.0
|
||||
pair_modify shift yes
|
||||
|
||||
neighbor 0.5 bin
|
||||
neigh_modify delay 0
|
||||
|
||||
change_box all triclinic
|
||||
|
||||
# these commads show the different methods that may be used to impose
|
||||
# a constant stress through isotropic or anisotropic coupling
|
||||
fix 1 all npt/uef temp 0.722 0.722 0.5 iso 10 10 5 erate 0.5 0.5 ext z
|
||||
#fix 1 all npt/uef temp 0.722 0.722 0.5 z 10 10 5 erate 0.5 0.5 ext xyz
|
||||
|
||||
fix 2 all momentum 100 linear 1 1 1
|
||||
|
||||
#dump 1 all atom 25 dump.lammpstrj
|
||||
|
||||
#dump 2 all cfg/uef 25 dump.*.cfg mass type xs ys zs
|
||||
|
||||
thermo_style custom step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
|
||||
thermo 50
|
||||
run 10000
|
||||
284
examples/USER/uef/npt_biaxial/log.22Sep17.npt_biaxial.g++.1
Normal file
284
examples/USER/uef/npt_biaxial/log.22Sep17.npt_biaxial.g++.1
Normal file
@ -0,0 +1,284 @@
|
||||
LAMMPS (22 Sep 2017)
|
||||
# biaxial NPT deformation of WCA fluid
|
||||
|
||||
units lj
|
||||
atom_style atomic
|
||||
|
||||
|
||||
pair_style lj/cut 1.122562
|
||||
read_data data.wca
|
||||
orthogonal box = (0 0 0) to (8.39798 8.39798 8.39798)
|
||||
1 by 1 by 1 MPI processor grid
|
||||
reading atoms ...
|
||||
500 atoms
|
||||
reading velocities ...
|
||||
500 velocities
|
||||
pair_coeff 1 1 1.0 1.0
|
||||
pair_modify shift yes
|
||||
|
||||
neighbor 0.5 bin
|
||||
neigh_modify delay 0
|
||||
|
||||
change_box all triclinic
|
||||
triclinic box = (0 0 0) to (8.39798 8.39798 8.39798) with tilt (0 0 0)
|
||||
|
||||
# these commads show the different methods that may be used to impose
|
||||
# a constant stress through isotropic or anisotropic coupling
|
||||
fix 1 all npt/uef temp 0.722 0.722 0.5 iso 10 10 5 erate 0.5 0.5 ext z
|
||||
#fix 1 all npt/uef temp 0.722 0.722 0.5 z 10 10 5 erate 0.5 0.5 ext xyz
|
||||
|
||||
fix 2 all momentum 100 linear 1 1 1
|
||||
|
||||
#dump 1 all atom 25 dump.lammpstrj
|
||||
|
||||
#dump 2 all cfg/uef 25 dump.*.cfg mass type xs ys zs
|
||||
|
||||
thermo_style custom step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
|
||||
thermo 50
|
||||
run 10000
|
||||
Neighbor list info ...
|
||||
update every 1 steps, delay 0 steps, check yes
|
||||
max neighbors/atom: 2000, page size: 100000
|
||||
master list distance cutoff = 1.62256
|
||||
ghost atom cutoff = 1.62256
|
||||
binsize = 0.811281, bins = 11 11 11
|
||||
1 neighbor lists, perpetual/occasional/extra = 1 0 0
|
||||
(1) pair lj/cut, perpetual
|
||||
attributes: half, newton on
|
||||
pair build: half/bin/newton/tri
|
||||
stencil: half/bin/3d/newton/tri
|
||||
bin: standard
|
||||
Per MPI rank memory allocation (min/avg/max) = 2.742 | 2.742 | 2.742 Mbytes
|
||||
Step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
0 6.3937851 7.0436438 6.4461087
|
||||
50 7.9020345 7.303289 14.733929
|
||||
100 8.3214325 8.385843 13.166247
|
||||
150 5.7168419 5.6083988 9.7224198
|
||||
200 3.8875348 4.0840705 7.9912973
|
||||
250 4.2933724 4.2410519 7.7200953
|
||||
300 4.5452314 4.2434949 8.6709832
|
||||
350 5.446489 5.6166962 11.300536
|
||||
400 6.0682558 5.7787878 12.106852
|
||||
450 5.3884296 5.5066688 9.9100012
|
||||
500 4.8046418 4.7115986 9.9769201
|
||||
550 4.9043933 4.6240966 9.319336
|
||||
600 4.6843238 4.9571246 9.5325487
|
||||
650 5.1952989 5.195648 9.6187047
|
||||
700 5.4163364 5.2938289 10.136655
|
||||
750 5.6723178 5.0670261 11.415221
|
||||
800 5.9379901 5.6402916 10.815209
|
||||
850 5.0695389 5.3021432 10.742859
|
||||
900 4.6498962 4.7111912 9.8453523
|
||||
950 4.4811641 5.170132 9.421765
|
||||
1000 4.7501483 4.2420772 9.3510803
|
||||
1050 4.5973379 5.2046799 9.8632975
|
||||
1100 4.7879562 4.9051316 9.8658626
|
||||
1150 5.0528771 5.5048545 10.110913
|
||||
1200 4.9926563 5.2482081 9.9610863
|
||||
1250 4.8476232 4.8498859 9.5752989
|
||||
1300 5.2932038 5.0594534 10.869545
|
||||
1350 5.4720421 5.0830442 10.69386
|
||||
1400 4.8175566 4.9863651 10.041898
|
||||
1450 4.6307145 4.5615459 9.231299
|
||||
1500 5.0296259 4.5433558 8.7180513
|
||||
1550 4.8708444 5.1601014 10.256791
|
||||
1600 5.5098593 5.316773 10.224386
|
||||
1650 5.5055907 5.0516814 10.272037
|
||||
1700 4.6965552 5.2436628 9.8272905
|
||||
1750 5.0212307 4.740601 8.9738802
|
||||
1800 4.7756351 5.199734 9.9929449
|
||||
1850 4.7107092 5.177203 10.580427
|
||||
1900 5.0935608 4.5286844 9.3036832
|
||||
1950 4.8035855 4.894359 9.4423987
|
||||
2000 4.5643937 4.7480477 9.4384251
|
||||
2050 4.9314701 4.990414 10.151815
|
||||
2100 5.9828824 5.8188026 12.262691
|
||||
2150 5.5337303 5.4598468 10.136112
|
||||
2200 4.892172 4.8699674 9.1629905
|
||||
2250 4.5680591 4.5740533 9.5633545
|
||||
2300 5.0023535 4.3948659 8.9645774
|
||||
2350 5.2533056 4.9803884 10.255653
|
||||
2400 5.3330196 5.3888322 10.021617
|
||||
2450 5.2095527 4.8765336 10.135381
|
||||
2500 5.221153 5.2974568 10.2493
|
||||
2550 5.385331 5.0801192 10.490479
|
||||
2600 5.3274988 5.0253548 10.147587
|
||||
2650 4.718677 5.2710337 9.7549521
|
||||
2700 4.5811521 4.6083971 8.1923164
|
||||
2750 4.4743752 4.3319737 8.8690805
|
||||
2800 5.0215013 5.2262961 9.8627954
|
||||
2850 6.1005107 5.5994942 11.170661
|
||||
2900 5.9673524 5.362124 10.401699
|
||||
2950 5.5400849 5.7523993 10.292536
|
||||
3000 4.9202636 5.1210431 10.125839
|
||||
3050 4.5286487 5.3586164 9.8634322
|
||||
3100 4.5719961 5.1615414 9.8027972
|
||||
3150 5.4091919 4.8581943 9.96938
|
||||
3200 5.522125 5.3564838 9.7638407
|
||||
3250 4.9019062 5.2514758 9.2993079
|
||||
3300 5.2375918 5.1439012 9.4313575
|
||||
3350 4.750235 4.8692016 10.54193
|
||||
3400 5.5793211 5.4184157 11.021389
|
||||
3450 4.9022614 5.3166498 9.4629659
|
||||
3500 4.6348617 4.4650979 9.0665548
|
||||
3550 4.7922405 4.8961269 10.255446
|
||||
3600 4.8914457 5.1158894 9.4736084
|
||||
3650 5.062771 4.6725475 10.263484
|
||||
3700 5.4842823 5.7793971 10.342915
|
||||
3750 5.3136012 5.063065 10.398307
|
||||
3800 4.9372149 4.9270414 9.5304748
|
||||
3850 5.2144752 5.1716455 9.7575725
|
||||
3900 5.0892665 5.1697057 9.918052
|
||||
3950 5.1124507 5.354702 9.791366
|
||||
4000 5.1255084 5.1143653 10.913101
|
||||
4050 5.1891698 4.9913681 9.6871728
|
||||
4100 4.7663368 4.2435014 8.3815668
|
||||
4150 4.8060033 4.3415868 9.6553386
|
||||
4200 4.8548303 4.8006768 9.5995801
|
||||
4250 5.0976366 5.2683175 10.386444
|
||||
4300 5.8921937 5.5134696 10.788143
|
||||
4350 5.8323871 5.5255869 11.199128
|
||||
4400 5.2464655 5.0005905 10.311055
|
||||
4450 4.9264849 5.2499854 10.26702
|
||||
4500 4.4431895 4.536981 8.7489096
|
||||
4550 4.5180914 4.2080277 8.6525529
|
||||
4600 5.1782188 5.1224059 10.683341
|
||||
4650 5.4156233 4.8714464 10.473939
|
||||
4700 5.3107669 5.224614 10.569391
|
||||
4750 4.9538022 5.2509475 10.288918
|
||||
4800 4.6976945 4.8107142 9.8299772
|
||||
4850 5.1227936 5.0737571 10.440452
|
||||
4900 4.7580514 4.6375995 9.1971008
|
||||
4950 5.0647601 4.6470735 9.583131
|
||||
5000 5.196231 5.7531491 10.409807
|
||||
5050 5.6691323 5.7163652 12.335701
|
||||
5100 5.3603738 5.4887106 10.961712
|
||||
5150 4.455028 4.6494465 9.8096968
|
||||
5200 4.7596912 4.4804896 9.3762885
|
||||
5250 5.3144927 5.0113772 9.553101
|
||||
5300 5.3445266 4.8262035 9.1220802
|
||||
5350 5.1540657 5.5982676 10.765178
|
||||
5400 5.1773418 5.2684381 10.452351
|
||||
5450 4.8946859 5.3283747 9.8015564
|
||||
5500 5.2009608 4.7183522 9.4558009
|
||||
5550 5.4158589 5.5005458 10.539505
|
||||
5600 4.7196831 5.4181991 9.6439249
|
||||
5650 4.8333571 4.8601728 8.9350189
|
||||
5700 5.4395698 4.9730096 10.669681
|
||||
5750 5.2947443 5.6973259 10.020539
|
||||
5800 5.4391993 5.5255143 10.264969
|
||||
5850 4.9921388 5.2643827 10.217028
|
||||
5900 5.0048643 4.7952641 8.9718929
|
||||
5950 5.1843818 4.5987295 9.6858944
|
||||
6000 5.0343993 4.946933 9.7436708
|
||||
6050 4.6202714 5.3502658 10.752915
|
||||
6100 5.6914422 5.3621964 10.281827
|
||||
6150 5.1928763 5.9652686 10.923881
|
||||
6200 5.0030409 5.2013891 10.056308
|
||||
6250 4.9699876 5.2363753 9.9964211
|
||||
6300 4.9129606 4.4558458 9.0419952
|
||||
6350 4.6623958 4.4078264 8.528649
|
||||
6400 4.9811441 5.1183207 10.261751
|
||||
6450 5.3644017 5.5153937 10.401295
|
||||
6500 5.6674981 5.7427566 11.928777
|
||||
6550 5.1622364 5.3212928 10.067198
|
||||
6600 4.5954278 5.1645397 10.16724
|
||||
6650 4.9192712 5.0413326 9.95656
|
||||
6700 4.6179845 4.5656214 9.3798952
|
||||
6750 4.7287495 4.5071148 8.7890116
|
||||
6800 4.8600442 4.8083512 10.245595
|
||||
6850 5.0515531 5.1609272 10.553855
|
||||
6900 5.1159742 5.1359869 10.594699
|
||||
6950 4.8908884 5.0592418 9.5698704
|
||||
7000 4.7654136 4.7530776 8.9439321
|
||||
7050 4.779293 4.7534957 9.7244349
|
||||
7100 5.2265831 5.6869073 10.32717
|
||||
7150 5.4019177 5.15174 10.457567
|
||||
7200 4.9817102 5.0596098 10.337574
|
||||
7250 5.1836654 5.6065238 10.723108
|
||||
7300 4.2916569 4.457143 8.5419099
|
||||
7350 4.3906104 4.5439294 9.0805415
|
||||
7400 4.998572 5.3386063 10.491418
|
||||
7450 5.1109022 5.0506801 10.636116
|
||||
7500 5.0248381 5.019932 10.217023
|
||||
7550 5.0109265 5.1438717 9.9032426
|
||||
7600 4.6628614 4.6204146 8.9459669
|
||||
7650 4.8930717 5.0650009 10.049331
|
||||
7700 4.9373454 5.6265835 10.210644
|
||||
7750 5.5001067 5.3133253 10.667995
|
||||
7800 5.0816102 5.0125753 10.591986
|
||||
7850 4.9638046 5.1969015 9.9728333
|
||||
7900 4.8438207 4.9217213 8.9978809
|
||||
7950 4.7318805 4.6248537 8.6806596
|
||||
8000 5.2808543 5.2892613 10.932535
|
||||
8050 5.9609722 5.87087 10.47602
|
||||
8100 5.2190231 5.6693244 11.244536
|
||||
8150 5.3481127 5.2849903 10.222845
|
||||
8200 4.7833053 4.7404657 9.2034474
|
||||
8250 4.5461994 4.510467 10.294452
|
||||
8300 4.6025175 4.8332817 8.7967546
|
||||
8350 5.0389897 5.6642908 10.243402
|
||||
8400 4.8414338 4.8925143 9.3653631
|
||||
8450 5.5087429 4.7830361 10.831666
|
||||
8500 5.2678146 5.1697789 9.9105782
|
||||
8550 5.1211843 4.9097801 9.4165956
|
||||
8600 5.8239149 5.0821022 10.803261
|
||||
8650 5.3620154 5.5831747 11.16202
|
||||
8700 5.1625813 4.8791404 10.537681
|
||||
8750 4.5622461 5.0157549 10.013227
|
||||
8800 4.4051517 5.0224553 9.6364273
|
||||
8850 4.1711629 4.635617 8.5470244
|
||||
8900 4.7049907 5.2458435 10.100728
|
||||
8950 4.8568883 5.2360772 9.2306469
|
||||
9000 5.0091899 5.2203574 10.718541
|
||||
9050 5.1037824 4.9022451 10.24271
|
||||
9100 5.0789015 4.9331454 9.173614
|
||||
9150 5.3865455 5.3427553 11.40199
|
||||
9200 5.5089482 5.9423232 10.976063
|
||||
9250 5.1353552 5.0650262 10.040607
|
||||
9300 4.6761948 4.9155175 9.6413722
|
||||
9350 4.4780834 4.3934708 8.7049819
|
||||
9400 4.2561799 4.7906324 9.046134
|
||||
9450 5.6162819 5.2881846 9.9040868
|
||||
9500 5.7554547 5.6111262 10.23849
|
||||
9550 5.4230462 5.5656045 10.908006
|
||||
9600 5.5045685 4.9818892 9.8929535
|
||||
9650 5.0541481 5.0183351 9.5226021
|
||||
9700 4.9712829 5.2395398 9.9996693
|
||||
9750 5.0960017 5.4419775 10.914719
|
||||
9800 5.0790688 5.6378474 10.00789
|
||||
9850 4.9661747 5.114502 9.4585052
|
||||
9900 5.0133498 4.7456254 9.4572653
|
||||
9950 5.3318846 4.6643122 10.096292
|
||||
10000 5.2227687 4.8924305 9.5894615
|
||||
Loop time of 4.78247 on 1 procs for 10000 steps with 500 atoms
|
||||
|
||||
Performance: 903298.340 tau/day, 2090.968 timesteps/s
|
||||
99.2% CPU use with 1 MPI tasks x no OpenMP threads
|
||||
|
||||
MPI task timing breakdown:
|
||||
Section | min time | avg time | max time |%varavg| %total
|
||||
---------------------------------------------------------------
|
||||
Pair | 0.73894 | 0.73894 | 0.73894 | 0.0 | 15.45
|
||||
Neigh | 2.9092 | 2.9092 | 2.9092 | 0.0 | 60.83
|
||||
Comm | 0.32306 | 0.32306 | 0.32306 | 0.0 | 6.76
|
||||
Output | 0.003392 | 0.003392 | 0.003392 | 0.0 | 0.07
|
||||
Modify | 0.6959 | 0.6959 | 0.6959 | 0.0 | 14.55
|
||||
Other | | 0.112 | | | 2.34
|
||||
|
||||
Nlocal: 500 ave 500 max 500 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
Nghost: 905 ave 905 max 905 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
Neighs: 3340 ave 3340 max 3340 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
|
||||
Total # of neighbors = 3340
|
||||
Ave neighs/atom = 6.68
|
||||
Neighbor list builds = 4580
|
||||
Dangerous builds = 0
|
||||
|
||||
Please see the log.cite file for references relevant to this simulation
|
||||
|
||||
Total wall time: 0:00:04
|
||||
284
examples/USER/uef/npt_biaxial/log.22Sep17.npt_biaxial.g++.4
Normal file
284
examples/USER/uef/npt_biaxial/log.22Sep17.npt_biaxial.g++.4
Normal file
@ -0,0 +1,284 @@
|
||||
LAMMPS (22 Sep 2017)
|
||||
# biaxial NPT deformation of WCA fluid
|
||||
|
||||
units lj
|
||||
atom_style atomic
|
||||
|
||||
|
||||
pair_style lj/cut 1.122562
|
||||
read_data data.wca
|
||||
orthogonal box = (0 0 0) to (8.39798 8.39798 8.39798)
|
||||
1 by 2 by 2 MPI processor grid
|
||||
reading atoms ...
|
||||
500 atoms
|
||||
reading velocities ...
|
||||
500 velocities
|
||||
pair_coeff 1 1 1.0 1.0
|
||||
pair_modify shift yes
|
||||
|
||||
neighbor 0.5 bin
|
||||
neigh_modify delay 0
|
||||
|
||||
change_box all triclinic
|
||||
triclinic box = (0 0 0) to (8.39798 8.39798 8.39798) with tilt (0 0 0)
|
||||
|
||||
# these commads show the different methods that may be used to impose
|
||||
# a constant stress through isotropic or anisotropic coupling
|
||||
fix 1 all npt/uef temp 0.722 0.722 0.5 iso 10 10 5 erate 0.5 0.5 ext z
|
||||
#fix 1 all npt/uef temp 0.722 0.722 0.5 z 10 10 5 erate 0.5 0.5 ext xyz
|
||||
|
||||
fix 2 all momentum 100 linear 1 1 1
|
||||
|
||||
#dump 1 all atom 25 dump.lammpstrj
|
||||
|
||||
#dump 2 all cfg/uef 25 dump.*.cfg mass type xs ys zs
|
||||
|
||||
thermo_style custom step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
|
||||
thermo 50
|
||||
run 10000
|
||||
Neighbor list info ...
|
||||
update every 1 steps, delay 0 steps, check yes
|
||||
max neighbors/atom: 2000, page size: 100000
|
||||
master list distance cutoff = 1.62256
|
||||
ghost atom cutoff = 1.62256
|
||||
binsize = 0.811281, bins = 11 11 11
|
||||
1 neighbor lists, perpetual/occasional/extra = 1 0 0
|
||||
(1) pair lj/cut, perpetual
|
||||
attributes: half, newton on
|
||||
pair build: half/bin/newton/tri
|
||||
stencil: half/bin/3d/newton/tri
|
||||
bin: standard
|
||||
Per MPI rank memory allocation (min/avg/max) = 2.71 | 2.71 | 2.71 Mbytes
|
||||
Step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
0 6.3937851 7.0436438 6.4461087
|
||||
50 7.9020345 7.303289 14.733929
|
||||
100 8.3214325 8.385843 13.166247
|
||||
150 5.7168419 5.6083988 9.7224198
|
||||
200 3.8875348 4.0840705 7.9912973
|
||||
250 4.2933724 4.2410519 7.7200953
|
||||
300 4.5452314 4.2434949 8.6709832
|
||||
350 5.446489 5.6166962 11.300536
|
||||
400 6.0682558 5.7787878 12.106852
|
||||
450 5.3884296 5.5066688 9.9100012
|
||||
500 4.8046418 4.7115986 9.9769201
|
||||
550 4.9043933 4.6240966 9.319336
|
||||
600 4.6843238 4.9571246 9.5325487
|
||||
650 5.1952989 5.195648 9.6187047
|
||||
700 5.4163364 5.2938289 10.136655
|
||||
750 5.6723178 5.0670261 11.415221
|
||||
800 5.9379901 5.6402916 10.815209
|
||||
850 5.0695389 5.3021432 10.742859
|
||||
900 4.6498961 4.7111912 9.8453524
|
||||
950 4.4811637 5.1701321 9.4217645
|
||||
1000 4.750149 4.2420768 9.3510788
|
||||
1050 4.5973376 5.2046787 9.8633025
|
||||
1100 4.7879517 4.9051339 9.8658578
|
||||
1150 5.0528775 5.5048671 10.110905
|
||||
1200 4.9926841 5.2482049 9.9610519
|
||||
1250 4.8475836 4.8499116 9.5753651
|
||||
1300 5.2930219 5.0593566 10.869335
|
||||
1350 5.4722342 5.0830411 10.693439
|
||||
1400 4.8165803 4.9851498 10.04213
|
||||
1450 4.6276458 4.5642988 9.2306141
|
||||
1500 5.0196773 4.5470773 8.7204145
|
||||
1550 4.878246 5.1583406 10.263895
|
||||
1600 5.4921049 5.3126759 10.274755
|
||||
1650 5.4363266 5.1708866 9.9880665
|
||||
1700 4.9186005 5.2841155 9.5911223
|
||||
1750 4.9105668 4.7112031 8.9221295
|
||||
1800 4.9833291 4.886821 9.6573393
|
||||
1850 5.0729703 4.8331712 10.094971
|
||||
1900 5.7220173 5.9330299 10.580261
|
||||
1950 5.3176022 5.7140521 11.11604
|
||||
2000 5.3247727 5.5172893 10.622834
|
||||
2050 5.2055235 4.8768078 9.9819356
|
||||
2100 4.4604981 4.1427844 7.8106592
|
||||
2150 4.4406592 4.8264893 9.6662695
|
||||
2200 5.2350836 5.1039145 10.36006
|
||||
2250 5.3777857 5.3274609 11.357157
|
||||
2300 5.4888334 5.9555482 10.76346
|
||||
2350 4.6122564 4.7356468 9.0833059
|
||||
2400 4.6670237 4.4895588 9.0619648
|
||||
2450 4.2201177 4.2558397 9.4898835
|
||||
2500 5.452448 5.4336384 10.50224
|
||||
2550 5.012581 5.316158 10.324517
|
||||
2600 5.0880279 5.1264772 10.085103
|
||||
2650 4.8885834 5.2368982 9.6002032
|
||||
2700 5.1549266 5.3419678 11.335447
|
||||
2750 5.497105 5.3643445 9.9990393
|
||||
2800 4.8826744 4.9875712 10.125435
|
||||
2850 4.8617121 5.1282348 9.5629401
|
||||
2900 4.7883446 4.6187804 9.0562496
|
||||
2950 4.7656266 5.1293592 10.693811
|
||||
3000 5.2676745 5.110172 9.3512146
|
||||
3050 4.5749222 4.8413907 10.311305
|
||||
3100 5.0794819 5.265009 9.1598496
|
||||
3150 5.2078869 5.1879882 10.412548
|
||||
3200 5.0187616 4.6226213 9.7266174
|
||||
3250 5.1541897 4.5157063 9.8355764
|
||||
3300 5.0721396 5.3545282 10.174356
|
||||
3350 5.3984495 5.3222207 10.008886
|
||||
3400 5.1263226 5.1189192 10.361534
|
||||
3450 5.1251845 4.8312752 9.6546597
|
||||
3500 5.1133696 5.2646289 10.320765
|
||||
3550 4.9884235 5.3861707 9.1944042
|
||||
3600 5.196909 5.203186 10.085965
|
||||
3650 5.4717592 5.2205442 10.251283
|
||||
3700 5.4429771 5.3027898 11.385714
|
||||
3750 5.5688484 5.5980199 10.558193
|
||||
3800 4.5239453 4.7021545 8.952588
|
||||
3850 4.6438079 4.6409958 9.3890154
|
||||
3900 5.1108473 4.8787691 10.665694
|
||||
3950 5.6398426 5.4386578 10.668189
|
||||
4000 5.063697 4.9663173 10.513266
|
||||
4050 4.8770847 4.4603573 9.8101845
|
||||
4100 4.3950768 4.3579384 8.3402845
|
||||
4150 4.3355402 5.0429352 10.323111
|
||||
4200 4.7688478 5.051487 9.0632339
|
||||
4250 4.9879366 5.3367146 10.409554
|
||||
4300 5.4578199 5.4889206 10.418789
|
||||
4350 5.6598068 5.4538572 10.842349
|
||||
4400 5.3705312 5.3796871 10.430547
|
||||
4450 4.663804 5.058851 9.2705923
|
||||
4500 4.3439039 4.3523422 8.1747925
|
||||
4550 4.5414802 4.3750772 9.2702452
|
||||
4600 4.9216199 5.2897069 10.747727
|
||||
4650 5.5154852 5.9628437 10.5168
|
||||
4700 5.45199 5.382787 10.654544
|
||||
4750 4.7525419 5.4701385 9.3189378
|
||||
4800 5.3696365 4.6134207 9.4455676
|
||||
4850 5.2444123 5.035993 9.4148435
|
||||
4900 5.6006507 4.8536828 10.283579
|
||||
4950 5.155711 4.978634 10.58973
|
||||
5000 5.0854607 4.9853307 9.2414296
|
||||
5050 5.1098462 4.7349164 9.8739001
|
||||
5100 5.1989395 5.0217416 9.8780949
|
||||
5150 5.612116 5.2165007 10.338464
|
||||
5200 5.0571356 5.3109846 10.685262
|
||||
5250 5.4832657 5.0371665 8.9420853
|
||||
5300 4.5312549 4.9629392 8.2478064
|
||||
5350 5.1617038 5.0533699 10.452218
|
||||
5400 5.7873394 5.6776926 11.926526
|
||||
5450 5.7002516 5.243239 10.940265
|
||||
5500 4.7896799 4.898544 10.163856
|
||||
5550 4.9155627 4.9567495 9.4445476
|
||||
5600 4.2447343 4.5045912 8.7732992
|
||||
5650 5.070197 4.7343938 9.9908239
|
||||
5700 4.9609446 5.0901934 10.812786
|
||||
5750 5.4001631 5.5552888 10.085896
|
||||
5800 5.4209837 4.7153245 9.6865245
|
||||
5850 4.9801041 5.180338 9.8930439
|
||||
5900 5.3423595 5.2341361 10.294159
|
||||
5950 5.683047 5.6830131 10.24313
|
||||
6000 5.0618789 5.4533644 10.713412
|
||||
6050 5.4034888 4.6341621 10.031976
|
||||
6100 5.1934299 4.7525347 9.1287151
|
||||
6150 5.0092398 4.806931 10.024305
|
||||
6200 5.3046516 5.3083532 9.6396223
|
||||
6250 5.2824046 4.7957557 10.305279
|
||||
6300 5.3007029 5.0071874 11.175322
|
||||
6350 5.1128883 4.990408 9.3439118
|
||||
6400 5.0543602 4.9971378 9.8259954
|
||||
6450 4.8843692 4.9116343 10.08132
|
||||
6500 4.5966453 4.8042861 9.160272
|
||||
6550 4.8510961 4.7096646 9.8009968
|
||||
6600 5.658307 5.2330511 10.739519
|
||||
6650 5.2374409 5.3241249 10.291779
|
||||
6700 4.9006975 5.0036186 9.9872029
|
||||
6750 5.2209104 5.1826025 9.5671875
|
||||
6800 4.5238727 4.3070529 8.6072303
|
||||
6850 4.2617247 4.7551571 9.7302077
|
||||
6900 5.6499354 4.8714257 10.723511
|
||||
6950 5.6881769 5.1800721 10.18134
|
||||
7000 5.588834 5.0104896 10.304105
|
||||
7050 4.9404045 4.8589121 9.7096741
|
||||
7100 5.2208179 4.9339808 9.7737491
|
||||
7150 5.4507842 5.046485 10.734783
|
||||
7200 4.9737171 5.5760486 9.1627431
|
||||
7250 4.5967409 4.750471 9.315832
|
||||
7300 5.5147308 5.3202861 10.542679
|
||||
7350 5.7730418 5.5363574 10.384376
|
||||
7400 4.9879586 5.2837443 9.4485798
|
||||
7450 5.1862431 4.8357842 10.017598
|
||||
7500 5.4528245 5.1864957 10.941774
|
||||
7550 5.4202434 5.1089468 10.128264
|
||||
7600 4.8063537 4.8723653 9.0364984
|
||||
7650 4.3144701 4.6148377 9.4939315
|
||||
7700 4.9033831 5.5327473 9.9054613
|
||||
7750 5.0693093 5.0768222 10.473081
|
||||
7800 5.0591805 5.6009473 10.006225
|
||||
7850 5.15269 5.468248 10.404619
|
||||
7900 5.1971759 5.0615117 9.9614488
|
||||
7950 4.9771238 4.886213 9.4730722
|
||||
8000 4.7731123 5.1111433 9.9550597
|
||||
8050 5.1655183 5.4432364 9.6649669
|
||||
8100 5.4597006 5.4026039 11.631184
|
||||
8150 5.3229643 5.4394219 9.9830611
|
||||
8200 4.5420712 4.9359646 9.0121988
|
||||
8250 4.5455108 4.3883996 10.304568
|
||||
8300 5.219721 5.1571958 10.305462
|
||||
8350 4.7291561 4.7391636 9.0768372
|
||||
8400 5.3262934 5.8221591 11.065466
|
||||
8450 4.7583026 5.2282086 10.291955
|
||||
8500 4.9174536 4.5701979 10.140444
|
||||
8550 4.459922 4.8293188 9.4438719
|
||||
8600 4.7962584 4.5811071 9.5158666
|
||||
8650 4.6097275 4.431952 9.4350505
|
||||
8700 5.166554 5.2000584 11.162202
|
||||
8750 5.2353596 5.1668944 10.829751
|
||||
8800 5.3150111 4.9983333 9.8402224
|
||||
8850 4.9571197 4.9196589 9.4880549
|
||||
8900 4.902223 4.7539187 10.002425
|
||||
8950 4.9531983 5.0517321 9.7838444
|
||||
9000 5.4003802 4.9900303 10.550808
|
||||
9050 4.9254643 5.0678701 10.24321
|
||||
9100 4.9902263 4.9056928 9.0326566
|
||||
9150 5.1003677 5.1555374 10.049056
|
||||
9200 5.2358131 5.5834504 10.354698
|
||||
9250 5.5781649 5.1188429 10.361369
|
||||
9300 5.2100192 5.0737267 10.074694
|
||||
9350 5.1462976 4.8010759 9.0279769
|
||||
9400 5.0177693 5.4890092 10.03612
|
||||
9450 5.2507957 5.541141 10.639854
|
||||
9500 5.836784 4.9498236 10.288015
|
||||
9550 5.4698183 5.8761209 10.979924
|
||||
9600 5.0534023 5.0043428 10.436057
|
||||
9650 4.3619773 4.5413125 9.047513
|
||||
9700 4.3777508 4.7902251 8.9501908
|
||||
9750 4.6851974 4.3152085 8.825764
|
||||
9800 5.7312665 4.850913 10.247637
|
||||
9850 6.2290614 5.5480801 10.934907
|
||||
9900 6.0316892 6.1098926 11.562223
|
||||
9950 5.6370814 5.3933342 11.148805
|
||||
10000 4.6429923 5.0853156 9.4267693
|
||||
Loop time of 3.06414 on 4 procs for 10000 steps with 500 atoms
|
||||
|
||||
Performance: 1409858.531 tau/day, 3263.561 timesteps/s
|
||||
99.4% CPU use with 4 MPI tasks x no OpenMP threads
|
||||
|
||||
MPI task timing breakdown:
|
||||
Section | min time | avg time | max time |%varavg| %total
|
||||
---------------------------------------------------------------
|
||||
Pair | 0.27954 | 0.30916 | 0.3866 | 8.1 | 10.09
|
||||
Neigh | 1.0657 | 1.129 | 1.2853 | 8.5 | 36.85
|
||||
Comm | 0.56334 | 0.80526 | 0.89816 | 15.6 | 26.28
|
||||
Output | 0.0032749 | 0.0041527 | 0.006541 | 2.1 | 0.14
|
||||
Modify | 0.6522 | 0.67016 | 0.69843 | 2.1 | 21.87
|
||||
Other | | 0.1464 | | | 4.78
|
||||
|
||||
Nlocal: 125 ave 130 max 122 min
|
||||
Histogram: 2 0 0 0 0 1 0 0 0 1
|
||||
Nghost: 459.75 ave 470 max 448 min
|
||||
Histogram: 1 0 0 0 0 1 1 0 0 1
|
||||
Neighs: 837 ave 867 max 811 min
|
||||
Histogram: 2 0 0 0 0 0 0 0 1 1
|
||||
|
||||
Total # of neighbors = 3348
|
||||
Ave neighs/atom = 6.696
|
||||
Neighbor list builds = 4589
|
||||
Dangerous builds = 0
|
||||
|
||||
Please see the log.cite file for references relevant to this simulation
|
||||
|
||||
Total wall time: 0:00:03
|
||||
1022
examples/USER/uef/nvt_uniaxial/data.wca
Normal file
1022
examples/USER/uef/nvt_uniaxial/data.wca
Normal file
File diff suppressed because it is too large
Load Diff
28
examples/USER/uef/nvt_uniaxial/in.nvt_uniaxial
Normal file
28
examples/USER/uef/nvt_uniaxial/in.nvt_uniaxial
Normal file
@ -0,0 +1,28 @@
|
||||
# uniaxial NVT deformation of WCA fluid
|
||||
|
||||
units lj
|
||||
atom_style atomic
|
||||
|
||||
|
||||
pair_style lj/cut 1.122562
|
||||
read_data data.wca
|
||||
pair_coeff 1 1 1.0 1.0
|
||||
pair_modify shift yes
|
||||
|
||||
neighbor 0.5 bin
|
||||
neigh_modify delay 0
|
||||
|
||||
change_box all triclinic
|
||||
|
||||
fix 1 all nvt/uef temp 0.722 0.722 0.5 erate -0.5 -0.5
|
||||
|
||||
fix 2 all momentum 100 linear 1 1 1
|
||||
|
||||
#dump 1 all atom 25 dump.lammpstrj
|
||||
|
||||
#dump 2 all cfg/uef 25 dump.*.cfg mass type xs ys zs
|
||||
|
||||
thermo_style custom step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
|
||||
thermo 50
|
||||
run 10000
|
||||
281
examples/USER/uef/nvt_uniaxial/log.22Sep17.npt_uniaxial.g++.1
Normal file
281
examples/USER/uef/nvt_uniaxial/log.22Sep17.npt_uniaxial.g++.1
Normal file
@ -0,0 +1,281 @@
|
||||
LAMMPS (22 Sep 2017)
|
||||
# uniaxial NVT deformation of WCA fluid
|
||||
|
||||
units lj
|
||||
atom_style atomic
|
||||
|
||||
|
||||
pair_style lj/cut 1.122562
|
||||
read_data data.wca
|
||||
orthogonal box = (0 0 0) to (8.39798 8.39798 8.39798)
|
||||
1 by 1 by 1 MPI processor grid
|
||||
reading atoms ...
|
||||
500 atoms
|
||||
reading velocities ...
|
||||
500 velocities
|
||||
pair_coeff 1 1 1.0 1.0
|
||||
pair_modify shift yes
|
||||
|
||||
neighbor 0.5 bin
|
||||
neigh_modify delay 0
|
||||
|
||||
change_box all triclinic
|
||||
triclinic box = (0 0 0) to (8.39798 8.39798 8.39798) with tilt (0 0 0)
|
||||
|
||||
fix 1 all nvt/uef temp 0.722 0.722 0.5 erate -0.5 -0.5
|
||||
|
||||
fix 2 all momentum 100 linear 1 1 1
|
||||
|
||||
#dump 1 all atom 25 dump.lammpstrj
|
||||
|
||||
#dump 2 all cfg/uef 25 dump.*.cfg mass type xs ys zs
|
||||
|
||||
thermo_style custom step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
|
||||
thermo 50
|
||||
run 10000
|
||||
Neighbor list info ...
|
||||
update every 1 steps, delay 0 steps, check yes
|
||||
max neighbors/atom: 2000, page size: 100000
|
||||
master list distance cutoff = 1.62256
|
||||
ghost atom cutoff = 1.62256
|
||||
binsize = 0.811281, bins = 11 11 11
|
||||
1 neighbor lists, perpetual/occasional/extra = 1 0 0
|
||||
(1) pair lj/cut, perpetual
|
||||
attributes: half, newton on
|
||||
pair build: half/bin/newton/tri
|
||||
stencil: half/bin/3d/newton/tri
|
||||
bin: standard
|
||||
Per MPI rank memory allocation (min/avg/max) = 2.742 | 2.742 | 2.742 Mbytes
|
||||
Step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
0 6.3937851 7.0436438 6.4461087
|
||||
50 10.369902 10.999889 6.5437384
|
||||
100 12.411384 13.145871 7.8892802
|
||||
150 12.88182 12.114068 7.5155182
|
||||
200 10.375571 10.979773 6.5624056
|
||||
250 10.158901 10.251273 5.4525068
|
||||
300 10.011314 9.8371382 5.4031393
|
||||
350 9.4890008 9.7992108 5.054963
|
||||
400 9.1715116 9.2647886 4.9600208
|
||||
450 9.8495961 9.3353483 4.986786
|
||||
500 9.5903809 8.576173 4.8252116
|
||||
550 8.8057509 9.7152078 4.9089022
|
||||
600 8.9743682 10.031277 4.8924072
|
||||
650 10.233852 9.3059527 5.1977058
|
||||
700 9.2361906 9.6439971 5.0923763
|
||||
750 10.449751 9.5911465 5.2548714
|
||||
800 9.8366857 9.9873885 4.9456956
|
||||
850 10.266351 9.9620134 5.0602713
|
||||
900 9.6325813 9.3910946 5.1466624
|
||||
950 9.1911919 9.1519006 4.8275345
|
||||
1000 9.2871712 9.8508501 5.2674099
|
||||
1050 9.7052035 10.291793 5.160359
|
||||
1100 9.9425552 10.361827 5.5520484
|
||||
1150 10.377563 10.26588 5.5997436
|
||||
1200 9.6368269 10.354413 5.4137873
|
||||
1250 9.2071055 9.1706334 4.9800411
|
||||
1300 9.7140516 8.7257777 4.8536943
|
||||
1350 10.145983 9.8247797 4.8921813
|
||||
1400 9.1998975 9.6935675 5.2269871
|
||||
1450 9.4838136 8.9035726 4.6867121
|
||||
1500 9.066389 9.811658 4.8875316
|
||||
1550 10.007083 9.8946118 5.204609
|
||||
1600 9.722137 9.9291588 5.326712
|
||||
1650 9.7424614 9.6995966 5.3577819
|
||||
1700 9.7185135 9.5922931 5.385856
|
||||
1750 9.6368481 8.1414594 5.0046079
|
||||
1800 9.8065782 9.6798036 4.92473
|
||||
1850 9.0463737 10.024205 5.6495695
|
||||
1900 9.0774484 10.527582 5.0363426
|
||||
1950 9.1691461 10.049703 4.8517514
|
||||
2000 9.7619602 10.049532 5.2088038
|
||||
2050 9.9203348 9.8056395 5.7699852
|
||||
2100 9.1455665 9.5864018 5.1410542
|
||||
2150 8.8501135 9.7917307 5.141824
|
||||
2200 9.1433685 9.8168548 5.1742441
|
||||
2250 9.3248968 10.454532 5.5456356
|
||||
2300 10.210505 9.361499 5.0935226
|
||||
2350 9.5041186 9.598986 5.2451005
|
||||
2400 8.968896 9.2811921 4.7799868
|
||||
2450 9.8201159 10.133892 4.9981817
|
||||
2500 9.4421084 9.7875766 5.4124809
|
||||
2550 9.6410643 9.9148309 5.0603288
|
||||
2600 9.3810676 9.5346849 5.3233836
|
||||
2650 9.2340337 8.7233538 5.0361302
|
||||
2700 9.405597 9.7848488 4.7537042
|
||||
2750 9.9447412 9.2366502 4.9736459
|
||||
2800 10.189443 9.452684 5.624362
|
||||
2850 9.6680124 9.0941543 5.0101469
|
||||
2900 8.4350889 9.1595018 4.7706866
|
||||
2950 9.0056117 9.3498593 4.7655963
|
||||
3000 9.4795416 9.7400747 5.4705868
|
||||
3050 9.6605264 9.4558374 5.1240166
|
||||
3100 9.681451 9.7966554 5.3174458
|
||||
3150 9.7887528 9.7371194 4.9349103
|
||||
3200 9.9720054 9.5592538 5.2306417
|
||||
3250 10.185023 9.9436037 5.4127312
|
||||
3300 9.9670676 9.6069307 5.026261
|
||||
3350 9.0090411 9.5975249 4.8881019
|
||||
3400 9.6535707 9.748683 4.9933873
|
||||
3450 9.7774977 9.7570511 4.8512619
|
||||
3500 9.7250053 10.36386 5.4057249
|
||||
3550 9.567788 9.773413 5.3936671
|
||||
3600 9.2094148 9.0211247 5.2328675
|
||||
3650 9.3512098 9.73407 4.576624
|
||||
3700 9.159074 9.2611591 5.2996448
|
||||
3750 9.4772798 8.9792211 5.1007079
|
||||
3800 9.3898357 9.1150609 5.2126508
|
||||
3850 9.2466312 9.1460651 4.867472
|
||||
3900 9.9316993 9.3061137 5.1219265
|
||||
3950 9.2550913 9.780254 5.371885
|
||||
4000 10.005356 9.7328972 5.2535791
|
||||
4050 9.5228242 9.5134113 5.4181393
|
||||
4100 9.576808 10.082425 5.1272491
|
||||
4150 10.20265 9.6667058 5.529118
|
||||
4200 9.9443671 10.09427 5.6163734
|
||||
4250 9.4956875 9.9462818 5.0820085
|
||||
4300 10.350144 10.003702 5.1621505
|
||||
4350 10.119486 9.8691507 5.4913457
|
||||
4400 9.4991604 10.516185 5.2512264
|
||||
4450 9.0812854 9.2835603 5.1695714
|
||||
4500 9.7538478 8.5298834 4.5954607
|
||||
4550 9.8920394 9.3581299 5.202587
|
||||
4600 9.2694921 9.3338291 5.201958
|
||||
4650 10.925818 9.5062049 4.8694531
|
||||
4700 11.020014 10.130224 5.582409
|
||||
4750 9.5005142 9.3571561 5.456739
|
||||
4800 9.6670147 9.6628702 5.0451252
|
||||
4850 9.134147 9.1308788 4.7950594
|
||||
4900 9.7466206 8.6171443 4.8716275
|
||||
4950 9.0397505 9.1996036 5.2010502
|
||||
5000 9.6898652 9.8914655 5.3959279
|
||||
5050 10.06771 9.7523891 5.3690408
|
||||
5100 10.000963 9.6279379 5.4077384
|
||||
5150 9.8686159 10.179702 5.0686824
|
||||
5200 9.8866247 9.725152 5.4350049
|
||||
5250 9.2068346 10.214424 5.3187713
|
||||
5300 9.713992 9.8069045 5.496359
|
||||
5350 9.423852 9.364085 5.2144606
|
||||
5400 9.4350241 9.5584633 5.0339001
|
||||
5450 10.555124 10.784922 5.1938072
|
||||
5500 9.4147344 10.33187 5.4360602
|
||||
5550 9.8514653 9.6575827 5.4959779
|
||||
5600 9.3138107 9.6592624 4.941387
|
||||
5650 9.1224809 8.7112257 5.0435936
|
||||
5700 8.8289158 10.749686 4.8916132
|
||||
5750 9.7200279 10.030606 5.2033161
|
||||
5800 9.8439873 9.6289015 5.5131934
|
||||
5850 9.6257294 9.4128988 4.9196038
|
||||
5900 9.7490214 9.5776313 5.0301815
|
||||
5950 9.1430855 10.108944 5.1406243
|
||||
6000 9.3358644 9.5633737 4.9787073
|
||||
6050 9.4432774 8.9464304 5.1466052
|
||||
6100 8.8878373 9.5048946 4.9190238
|
||||
6150 9.6451898 9.2419823 5.0159841
|
||||
6200 9.5042173 8.9414307 5.2634247
|
||||
6250 9.0896505 9.7230651 5.3340322
|
||||
6300 8.8100599 8.8781352 5.4117914
|
||||
6350 9.3104601 9.0233294 5.3136432
|
||||
6400 9.368101 9.6387362 4.7833216
|
||||
6450 10.334343 9.8384149 5.3606204
|
||||
6500 9.8824036 10.022627 6.0857086
|
||||
6550 9.7034443 10.026765 5.17604
|
||||
6600 9.3757845 9.899268 5.2301359
|
||||
6650 10.540821 10.4343 5.5287065
|
||||
6700 9.6317649 9.8923579 5.6045768
|
||||
6750 9.5982772 10.07897 5.1221451
|
||||
6800 10.239883 10.189967 5.3167447
|
||||
6850 10.017271 9.7680902 5.229621
|
||||
6900 9.6200416 10.129301 5.1998759
|
||||
6950 9.0361417 8.923798 5.1652612
|
||||
7000 9.3153521 9.063054 4.6860773
|
||||
7050 8.6434091 9.0363436 4.7811975
|
||||
7100 9.4955395 9.3830541 5.022538
|
||||
7150 9.3392402 9.1847119 5.1544622
|
||||
7200 9.4676321 9.8370036 4.8854442
|
||||
7250 9.5115882 10.350324 4.9780525
|
||||
7300 9.6025583 9.6247917 5.473794
|
||||
7350 9.8919524 10.049446 4.9816931
|
||||
7400 9.6814319 9.9410894 5.265078
|
||||
7450 9.4130955 10.191436 5.2531256
|
||||
7500 9.8114668 8.8461635 5.0562894
|
||||
7550 10.321567 9.4730124 5.2043655
|
||||
7600 9.5059024 9.8330367 5.0749721
|
||||
7650 10.067084 10.606423 5.5598818
|
||||
7700 10.896159 10.084281 5.5159718
|
||||
7750 9.754306 10.162301 5.2475876
|
||||
7800 9.7278145 9.801009 5.0685504
|
||||
7850 9.8639905 10.323104 5.2458864
|
||||
7900 9.7246799 9.1377357 5.1841319
|
||||
7950 10.381792 9.6977533 5.0977386
|
||||
8000 9.7265224 9.2999829 4.7925571
|
||||
8050 9.5203288 9.3144956 4.7539211
|
||||
8100 9.6900973 9.6614063 5.304712
|
||||
8150 9.713677 9.9154149 4.8178575
|
||||
8200 9.4733597 9.8948632 5.036112
|
||||
8250 9.7783036 9.9554334 5.3355682
|
||||
8300 9.4034783 9.9156801 5.5539279
|
||||
8350 9.2984025 9.2013949 5.0753991
|
||||
8400 9.9299078 9.6166801 5.115776
|
||||
8450 9.526737 9.3312125 4.7800587
|
||||
8500 9.581956 10.065906 4.9756092
|
||||
8550 9.2767953 9.326006 5.3024978
|
||||
8600 9.799968 8.4660845 5.0199109
|
||||
8650 9.8985354 10.127852 4.9098064
|
||||
8700 8.7952691 10.521133 5.6840528
|
||||
8750 9.8299997 9.5588553 5.3085734
|
||||
8800 9.0811776 9.5704532 5.1684993
|
||||
8850 9.8303571 9.7618932 5.1251259
|
||||
8900 9.9238794 9.9654863 5.2999683
|
||||
8950 10.851304 9.9682289 5.4133763
|
||||
9000 9.5523794 9.1890766 5.1469144
|
||||
9050 9.7461948 8.9611236 4.9490826
|
||||
9100 10.138917 9.6757567 5.0473544
|
||||
9150 9.4869835 9.4786575 5.0142464
|
||||
9200 10.263518 10.079135 5.1493398
|
||||
9250 9.8691684 9.0908275 5.3221203
|
||||
9300 9.8586707 9.4177643 5.1525265
|
||||
9350 9.3375816 9.9167208 5.4846207
|
||||
9400 9.5603903 9.4813199 4.6237495
|
||||
9450 10.30892 9.5217736 5.6163214
|
||||
9500 9.327949 9.9831649 4.8923915
|
||||
9550 9.8421656 9.3202702 5.3352046
|
||||
9600 8.8543704 9.4556702 4.6430041
|
||||
9650 10.103399 9.2161072 4.8658062
|
||||
9700 9.507811 9.9647378 5.227369
|
||||
9750 9.4988096 8.9942893 5.2491418
|
||||
9800 9.8007958 9.234452 5.1740203
|
||||
9850 9.6029685 10.076042 5.023107
|
||||
9900 9.4035691 10.13782 5.2775777
|
||||
9950 9.6517135 10.355994 5.3035779
|
||||
10000 9.9157616 9.7741952 5.5269431
|
||||
Loop time of 3.47119 on 1 procs for 10000 steps with 500 atoms
|
||||
|
||||
Performance: 1244529.812 tau/day, 2880.856 timesteps/s
|
||||
99.5% CPU use with 1 MPI tasks x no OpenMP threads
|
||||
|
||||
MPI task timing breakdown:
|
||||
Section | min time | avg time | max time |%varavg| %total
|
||||
---------------------------------------------------------------
|
||||
Pair | 0.53077 | 0.53077 | 0.53077 | 0.0 | 15.29
|
||||
Neigh | 2.219 | 2.219 | 2.219 | 0.0 | 63.93
|
||||
Comm | 0.21866 | 0.21866 | 0.21866 | 0.0 | 6.30
|
||||
Output | 0.0027909 | 0.0027909 | 0.0027909 | 0.0 | 0.08
|
||||
Modify | 0.41956 | 0.41956 | 0.41956 | 0.0 | 12.09
|
||||
Other | | 0.08042 | | | 2.32
|
||||
|
||||
Nlocal: 500 ave 500 max 500 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
Nghost: 908 ave 908 max 908 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
Neighs: 3489 ave 3489 max 3489 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
|
||||
Total # of neighbors = 3489
|
||||
Ave neighs/atom = 6.978
|
||||
Neighbor list builds = 4567
|
||||
Dangerous builds = 0
|
||||
|
||||
Please see the log.cite file for references relevant to this simulation
|
||||
|
||||
Total wall time: 0:00:03
|
||||
281
examples/USER/uef/nvt_uniaxial/log.22Sep17.npt_uniaxial.g++.4
Normal file
281
examples/USER/uef/nvt_uniaxial/log.22Sep17.npt_uniaxial.g++.4
Normal file
@ -0,0 +1,281 @@
|
||||
LAMMPS (22 Sep 2017)
|
||||
# uniaxial NVT deformation of WCA fluid
|
||||
|
||||
units lj
|
||||
atom_style atomic
|
||||
|
||||
|
||||
pair_style lj/cut 1.122562
|
||||
read_data data.wca
|
||||
orthogonal box = (0 0 0) to (8.39798 8.39798 8.39798)
|
||||
1 by 2 by 2 MPI processor grid
|
||||
reading atoms ...
|
||||
500 atoms
|
||||
reading velocities ...
|
||||
500 velocities
|
||||
pair_coeff 1 1 1.0 1.0
|
||||
pair_modify shift yes
|
||||
|
||||
neighbor 0.5 bin
|
||||
neigh_modify delay 0
|
||||
|
||||
change_box all triclinic
|
||||
triclinic box = (0 0 0) to (8.39798 8.39798 8.39798) with tilt (0 0 0)
|
||||
|
||||
fix 1 all nvt/uef temp 0.722 0.722 0.5 erate -0.5 -0.5
|
||||
|
||||
fix 2 all momentum 100 linear 1 1 1
|
||||
|
||||
#dump 1 all atom 25 dump.lammpstrj
|
||||
|
||||
#dump 2 all cfg/uef 25 dump.*.cfg mass type xs ys zs
|
||||
|
||||
thermo_style custom step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
|
||||
thermo 50
|
||||
run 10000
|
||||
Neighbor list info ...
|
||||
update every 1 steps, delay 0 steps, check yes
|
||||
max neighbors/atom: 2000, page size: 100000
|
||||
master list distance cutoff = 1.62256
|
||||
ghost atom cutoff = 1.62256
|
||||
binsize = 0.811281, bins = 11 11 11
|
||||
1 neighbor lists, perpetual/occasional/extra = 1 0 0
|
||||
(1) pair lj/cut, perpetual
|
||||
attributes: half, newton on
|
||||
pair build: half/bin/newton/tri
|
||||
stencil: half/bin/3d/newton/tri
|
||||
bin: standard
|
||||
Per MPI rank memory allocation (min/avg/max) = 2.71 | 2.71 | 2.71 Mbytes
|
||||
Step c_1_press[1] c_1_press[2] c_1_press[3]
|
||||
0 6.3937851 7.0436438 6.4461087
|
||||
50 10.369902 10.999889 6.5437384
|
||||
100 12.411384 13.145871 7.8892802
|
||||
150 12.88182 12.114068 7.5155182
|
||||
200 10.375571 10.979773 6.5624056
|
||||
250 10.158901 10.251273 5.4525068
|
||||
300 10.011314 9.8371382 5.4031393
|
||||
350 9.4890008 9.7992108 5.054963
|
||||
400 9.1715116 9.2647886 4.9600208
|
||||
450 9.8495961 9.3353483 4.986786
|
||||
500 9.5903809 8.576173 4.8252116
|
||||
550 8.8057509 9.7152078 4.9089022
|
||||
600 8.9743682 10.031277 4.8924072
|
||||
650 10.233852 9.3059527 5.1977058
|
||||
700 9.2361906 9.6439971 5.0923762
|
||||
750 10.449751 9.5911465 5.2548714
|
||||
800 9.8366858 9.9873885 4.9456956
|
||||
850 10.266351 9.9620134 5.0602713
|
||||
900 9.6325811 9.3910947 5.1466627
|
||||
950 9.1911906 9.1518996 4.8275336
|
||||
1000 9.2871726 9.8508561 5.2674079
|
||||
1050 9.7052019 10.291782 5.1603609
|
||||
1100 9.9425827 10.361848 5.5520588
|
||||
1150 10.37757 10.265891 5.5997246
|
||||
1200 9.6367069 10.354451 5.4138749
|
||||
1250 9.2075006 9.1705136 4.9806417
|
||||
1300 9.7126123 8.7247301 4.8528533
|
||||
1350 10.146034 9.8319205 4.8958266
|
||||
1400 9.19638 9.7052094 5.2332913
|
||||
1450 9.4667073 8.9014702 4.7073387
|
||||
1500 9.0949401 9.8381729 4.8904182
|
||||
1550 10.068959 9.8096692 5.1613658
|
||||
1600 10.131443 9.7086921 5.0913546
|
||||
1650 9.9557572 9.9183989 5.2635245
|
||||
1700 9.3820418 8.6135867 4.9395498
|
||||
1750 8.9486365 9.3433293 5.2674469
|
||||
1800 10.352322 10.354432 5.3335153
|
||||
1850 8.9629784 9.6172224 5.0824828
|
||||
1900 8.5429652 9.749815 5.1577895
|
||||
1950 8.9828002 8.9631646 5.0948426
|
||||
2000 9.4171989 9.7122117 5.1386251
|
||||
2050 9.767299 9.6146113 4.8965319
|
||||
2100 10.049399 10.02243 5.3856622
|
||||
2150 9.2639644 10.09669 5.3020322
|
||||
2200 9.1916162 10.841127 5.0271109
|
||||
2250 9.9859708 8.519568 5.24568
|
||||
2300 9.9462279 9.5850706 5.0632906
|
||||
2350 9.8043011 9.7605126 5.3148323
|
||||
2400 9.089818 9.6474522 5.0012486
|
||||
2450 9.655874 9.3255636 4.8291262
|
||||
2500 9.0259445 9.3074827 5.1593031
|
||||
2550 9.1610478 9.7356492 5.0741161
|
||||
2600 9.3153881 9.3936949 5.4830363
|
||||
2650 9.6212696 9.598252 4.7172875
|
||||
2700 9.318906 9.665656 5.1917527
|
||||
2750 9.6613838 9.7106796 5.1843146
|
||||
2800 10.231844 9.9407641 4.6940455
|
||||
2850 10.008093 9.2781076 5.2624957
|
||||
2900 10.464168 9.808666 5.1457978
|
||||
2950 9.7080601 9.6972304 5.406456
|
||||
3000 9.6851209 10.050737 5.1198394
|
||||
3050 8.9093581 9.213428 5.233108
|
||||
3100 8.8098921 9.6607476 4.9625331
|
||||
3150 8.6608386 9.7503441 5.0737533
|
||||
3200 9.7004403 9.6748778 5.0734462
|
||||
3250 10.077054 10.318711 5.3233841
|
||||
3300 10.63887 9.4901739 5.516542
|
||||
3350 10.232909 9.3407073 5.4989967
|
||||
3400 9.8854134 9.9542625 5.1982468
|
||||
3450 9.9388203 9.3394716 4.8135833
|
||||
3500 9.2514026 9.5857527 5.4605449
|
||||
3550 8.3940282 9.0549836 4.6951548
|
||||
3600 8.7696659 9.6262816 4.6302087
|
||||
3650 9.0695173 8.8520895 5.0814723
|
||||
3700 9.4700744 9.2777557 4.6825004
|
||||
3750 9.0221844 9.3407486 5.3984156
|
||||
3800 9.1478369 9.1420043 5.3024474
|
||||
3850 9.2652818 9.700053 5.3355364
|
||||
3900 10.248456 9.2480211 5.2375956
|
||||
3950 9.8259922 9.6130415 5.5978761
|
||||
4000 9.404877 9.7931698 5.3767927
|
||||
4050 10.131713 9.7047295 5.2964594
|
||||
4100 9.8128638 9.4253237 5.5308166
|
||||
4150 10.25183 9.1333595 5.1957555
|
||||
4200 10.449736 9.295762 4.7863033
|
||||
4250 9.7304858 9.4482515 5.356439
|
||||
4300 9.2773777 9.4110855 4.9879246
|
||||
4350 8.9325082 9.3429549 5.0410132
|
||||
4400 9.2603855 9.7905381 4.7436126
|
||||
4450 8.7600443 9.9160722 5.196316
|
||||
4500 9.0824514 10.036035 4.873051
|
||||
4550 9.3884333 9.6644343 5.1154951
|
||||
4600 9.8181676 9.6385118 5.3639835
|
||||
4650 9.1574799 9.5386974 4.7487094
|
||||
4700 10.275911 9.5383553 5.1084297
|
||||
4750 9.2476854 9.3353591 5.4773008
|
||||
4800 9.3422407 9.1931821 4.9210291
|
||||
4850 9.9033126 9.6443642 5.1334553
|
||||
4900 9.8469467 9.6836455 5.5101146
|
||||
4950 10.211229 9.760253 5.0151668
|
||||
5000 9.3256273 10.109873 5.4129479
|
||||
5050 9.9704879 9.8504809 5.2191163
|
||||
5100 9.4820718 9.711504 4.9633504
|
||||
5150 9.2952171 9.6955742 5.3214246
|
||||
5200 9.8400791 9.4763906 4.9149518
|
||||
5250 10.080112 9.5634741 5.1701348
|
||||
5300 9.4221014 9.8149742 5.3605931
|
||||
5350 8.8228402 10.121343 5.3192212
|
||||
5400 9.7318719 10.508763 5.0044083
|
||||
5450 8.7429847 9.6583774 5.033313
|
||||
5500 9.9243256 10.239521 5.2956506
|
||||
5550 9.6088558 9.8555986 5.359257
|
||||
5600 10.904352 10.72033 5.7874034
|
||||
5650 9.4774477 9.7246962 5.596949
|
||||
5700 9.8527139 9.2188245 4.6518586
|
||||
5750 9.7932375 9.2763721 4.6663307
|
||||
5800 9.5630073 9.1704583 5.4784197
|
||||
5850 10.159996 9.5206168 5.0012706
|
||||
5900 9.1667978 9.6891715 5.1959301
|
||||
5950 9.2194131 9.2170699 5.1653264
|
||||
6000 9.5917124 9.2038051 5.1004966
|
||||
6050 9.4141124 9.9166471 5.0535712
|
||||
6100 10.231166 9.7746591 5.2399634
|
||||
6150 9.6054192 9.34446 5.0711646
|
||||
6200 9.9279645 8.9546561 5.4698039
|
||||
6250 9.1581437 8.95439 4.9336111
|
||||
6300 8.9257232 9.0665473 5.188718
|
||||
6350 9.7685695 9.0822789 4.8454457
|
||||
6400 9.5140226 9.5073414 5.3420644
|
||||
6450 9.5379198 9.3316002 5.1940338
|
||||
6500 10.216584 9.8982859 5.2790157
|
||||
6550 10.221286 10.913585 5.4521496
|
||||
6600 10.925166 10.088055 5.2079758
|
||||
6650 10.012501 9.3870455 5.1740433
|
||||
6700 10.176139 9.4869588 5.3783369
|
||||
6750 9.2819415 8.1431975 4.8256476
|
||||
6800 9.790723 9.7932943 5.3602927
|
||||
6850 9.1327199 9.8354267 5.1356642
|
||||
6900 9.240909 9.3854618 5.1767268
|
||||
6950 10.652806 9.4355671 5.1314851
|
||||
7000 9.1676917 10.055501 5.2824651
|
||||
7050 9.6486488 9.1063089 4.9596976
|
||||
7100 9.3733678 9.5875095 5.0882285
|
||||
7150 9.282243 9.354934 4.6752644
|
||||
7200 8.7370777 10.1445 5.0217245
|
||||
7250 9.563419 9.929442 5.1720916
|
||||
7300 9.9128346 9.2426863 5.3595023
|
||||
7350 9.9319949 9.9479477 5.4844629
|
||||
7400 9.4671739 9.2380048 5.3274217
|
||||
7450 9.7531193 9.4075641 4.9426738
|
||||
7500 10.331422 9.4704554 5.3036636
|
||||
7550 8.7724434 9.663179 4.9492563
|
||||
7600 10.222748 9.1674852 4.8135992
|
||||
7650 9.5288583 8.5866929 5.0342955
|
||||
7700 8.8832847 8.8035037 4.6184473
|
||||
7750 9.4847076 9.7969656 5.0035855
|
||||
7800 9.2782486 9.1283641 5.3151299
|
||||
7850 8.962989 8.7729872 5.0374817
|
||||
7900 9.6336556 9.4975426 4.9832641
|
||||
7950 10.125279 10.044247 5.2740873
|
||||
8000 9.9705097 9.6925988 5.2344188
|
||||
8050 9.7836318 9.5794261 4.9174534
|
||||
8100 10.79561 10.04809 5.1604639
|
||||
8150 9.7261726 10.372878 5.3083843
|
||||
8200 10.338569 10.08189 5.6555172
|
||||
8250 9.7833315 9.9834209 5.3790753
|
||||
8300 9.5857944 9.6052323 5.0689289
|
||||
8350 9.6706688 9.6189242 4.7725678
|
||||
8400 10.570423 9.5786093 5.1894242
|
||||
8450 9.6514599 10.168359 5.0733592
|
||||
8500 10.273682 9.9179284 5.4014563
|
||||
8550 9.3120287 10.146837 4.9895115
|
||||
8600 9.511943 9.644112 5.462624
|
||||
8650 10.380674 9.1117114 5.156727
|
||||
8700 10.068596 8.7687113 5.1440814
|
||||
8750 9.2484971 9.2477678 4.9318794
|
||||
8800 9.7298469 8.9480303 5.1151321
|
||||
8850 9.7299502 10.415138 4.7902908
|
||||
8900 10.966912 10.732962 5.4793574
|
||||
8950 10.328384 9.9501313 5.6238396
|
||||
9000 9.7385041 9.8319224 5.1926497
|
||||
9050 9.7971055 9.5740203 5.1111302
|
||||
9100 9.7789727 9.9281901 5.1786549
|
||||
9150 9.9306964 9.3360599 4.9524547
|
||||
9200 9.8798841 10.240752 5.1691344
|
||||
9250 10.185445 9.4934917 4.9188964
|
||||
9300 8.9184663 8.9349408 4.8079511
|
||||
9350 9.6552187 9.9846949 4.9619969
|
||||
9400 10.304306 9.2298208 5.2822855
|
||||
9450 9.8379613 10.041703 5.4186514
|
||||
9500 10.221443 9.5342818 4.8929802
|
||||
9550 9.9723047 10.072856 5.4169676
|
||||
9600 9.3923879 9.2984387 5.4452785
|
||||
9650 8.9072589 9.7482374 4.7835208
|
||||
9700 9.8370121 10.205922 5.0385145
|
||||
9750 9.4274542 9.4653248 5.6340681
|
||||
9800 9.7668106 9.3265705 5.3154126
|
||||
9850 10.422549 10.362922 5.361592
|
||||
9900 9.6264407 9.9790162 5.6381052
|
||||
9950 10.35018 9.8853593 5.2639184
|
||||
10000 9.6190853 9.7903758 5.1583115
|
||||
Loop time of 1.32471 on 4 procs for 10000 steps with 500 atoms
|
||||
|
||||
Performance: 3261084.476 tau/day, 7548.807 timesteps/s
|
||||
98.4% CPU use with 4 MPI tasks x no OpenMP threads
|
||||
|
||||
MPI task timing breakdown:
|
||||
Section | min time | avg time | max time |%varavg| %total
|
||||
---------------------------------------------------------------
|
||||
Pair | 0.1301 | 0.13489 | 0.13886 | 0.9 | 10.18
|
||||
Neigh | 0.5654 | 0.575 | 0.5883 | 1.3 | 43.41
|
||||
Comm | 0.35135 | 0.37147 | 0.38856 | 2.4 | 28.04
|
||||
Output | 0.0028336 | 0.0032777 | 0.004292 | 1.0 | 0.25
|
||||
Modify | 0.16328 | 0.16711 | 0.17231 | 0.8 | 12.61
|
||||
Other | | 0.07297 | | | 5.51
|
||||
|
||||
Nlocal: 125 ave 128 max 122 min
|
||||
Histogram: 1 1 0 0 0 0 0 0 1 1
|
||||
Nghost: 477 ave 489 max 469 min
|
||||
Histogram: 2 0 0 0 0 1 0 0 0 1
|
||||
Neighs: 861 ave 877 max 848 min
|
||||
Histogram: 1 0 1 0 0 1 0 0 0 1
|
||||
|
||||
Total # of neighbors = 3444
|
||||
Ave neighs/atom = 6.888
|
||||
Neighbor list builds = 4560
|
||||
Dangerous builds = 0
|
||||
|
||||
Please see the log.cite file for references relevant to this simulation
|
||||
|
||||
Total wall time: 0:00:01
|
||||
9
src/.gitignore
vendored
9
src/.gitignore
vendored
@ -103,6 +103,15 @@
|
||||
/sna.cpp
|
||||
/sna.h
|
||||
|
||||
/uef_*.cpp
|
||||
/uef_*.h
|
||||
/compute_*uef.cpp
|
||||
/compute_*uef.h
|
||||
/dump_*uef.cpp
|
||||
/dump_*uef.h
|
||||
/fix_*uef.cpp
|
||||
/fix_*uef.h
|
||||
|
||||
/atom_vec_wavepacket.cpp
|
||||
/atom_vec_wavepacket.h
|
||||
/fix_nve_awpmd.cpp
|
||||
|
||||
@ -63,7 +63,7 @@ PACKUSER = user-atc user-awpmd user-cgdna user-cgsdk user-colvars \
|
||||
user-mgpt user-misc user-molfile \
|
||||
user-netcdf user-omp user-phonon user-qmmm user-qtb \
|
||||
user-quip user-reaxc user-smd user-smtbq user-sph user-tally \
|
||||
user-vtk
|
||||
user-uef user-vtk
|
||||
|
||||
PACKLIB = compress gpu kim kokkos latte meam mpiio mscg poems \
|
||||
python reax voronoi \
|
||||
|
||||
47
src/USER-UEF/README
Normal file
47
src/USER-UEF/README
Normal file
@ -0,0 +1,47 @@
|
||||
USER-UEF is a LAMMPS package for non-equilibrium molecular dynamics
|
||||
(NEMD) under diagonal flow fields, including uniaxial and biaxial
|
||||
flow. With this package, simulations under extensional flow may be
|
||||
carried out for an indefinite amount of time. It is an implementation
|
||||
of the boundary conditions developed by Matthew Dobson, and also uses
|
||||
numerical lattice reduction as was proposed by Thomas Hunt. The
|
||||
lattice reduction algorithm is from Igor Semaev. The package is
|
||||
intended for simulations of homogeneous flows, and integrates the
|
||||
SLLOD equations of motion.
|
||||
|
||||
--
|
||||
|
||||
This package was created by David Nicholson (davidanich@gmail.com)
|
||||
at the Massachusetts Institute of Technology.
|
||||
|
||||
--
|
||||
|
||||
The following commands are contained in this package:
|
||||
|
||||
fix npt/uef and fix nvt/uef:
|
||||
These commands perform time-integration of the SLLOD equations of
|
||||
motion under constant temperature/pressure with the proper
|
||||
boundary conditions for extensional flow fields.
|
||||
|
||||
compute pressure/uef and compute temp/uef:
|
||||
These commands allow for the evaluation of pressure/ke tensors in
|
||||
the reference frame corresponding to the applied flow field.
|
||||
|
||||
dump cfg/uef:
|
||||
This command dumps coordinates in the reference frame corresponding
|
||||
to the applied flow field.
|
||||
|
||||
For more information, visit the documentation page for fix nvt/uef
|
||||
and examine the example scripts in doc/USER/uef/.
|
||||
|
||||
--
|
||||
|
||||
References:
|
||||
|
||||
I. Semaev, Cryptography and Lattices, 181 (2001).
|
||||
|
||||
M. Dobson, J. Chem. Phys., 141, 184103 (2014)
|
||||
|
||||
T.A. Hunt, Mol. Simul., 42, 347 (2016).
|
||||
|
||||
D.A. Nicholson, G.C. Rutledge, J. Chem. Phys.,145, 244903 (2016)
|
||||
|
||||
195
src/USER-UEF/compute_pressure_uef.cpp
Normal file
195
src/USER-UEF/compute_pressure_uef.cpp
Normal file
@ -0,0 +1,195 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include <mpi.h>
|
||||
#include <string.h>
|
||||
#include <stdlib.h>
|
||||
#include "compute_pressure_uef.h"
|
||||
#include "fix_nh_uef.h"
|
||||
#include "update.h"
|
||||
#include "domain.h"
|
||||
#include "modify.h"
|
||||
#include "fix.h"
|
||||
#include "force.h"
|
||||
#include "pair.h"
|
||||
#include "bond.h"
|
||||
#include "angle.h"
|
||||
#include "dihedral.h"
|
||||
#include "improper.h"
|
||||
#include "kspace.h"
|
||||
#include "error.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Default values for the ext flags
|
||||
* ----------------------------------------------------------------------*/
|
||||
ComputePressureUef::ComputePressureUef(LAMMPS *lmp, int narg, char **arg) :
|
||||
ComputePressure(lmp, narg, arg)
|
||||
{
|
||||
ext_flags[0] = true;
|
||||
ext_flags[1] = true;
|
||||
ext_flags[2] = true;
|
||||
in_fix=false;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Check for the uef fix
|
||||
* ----------------------------------------------------------------------*/
|
||||
void ComputePressureUef::init()
|
||||
{
|
||||
ComputePressure::init();
|
||||
// check to make sure the other uef fix is on
|
||||
// borrowed from Pieter's nvt/sllod code
|
||||
int i=0;
|
||||
for (i=0; i<modify->nfix; i++) {
|
||||
if (strcmp(modify->fix[i]->style,"nvt/uef")==0)
|
||||
break;
|
||||
if (strcmp(modify->fix[i]->style,"npt/uef")==0)
|
||||
break;
|
||||
}
|
||||
if (i==modify->nfix)
|
||||
error->all(FLERR,"Can't use compute pressure/uef without defining a fix nvt/npt/uef");
|
||||
ifix_uef=i;
|
||||
((FixNHUef*) modify->fix[ifix_uef])->get_ext_flags(ext_flags);
|
||||
|
||||
if (strcmp(temperature->style,"temp/uef") != 0)
|
||||
error->warning(FLERR,"The temperature used in compute pressure/ued is not of style temp/uef");
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Compute pressure in the directions i corresponding to ext_flag[i]=true
|
||||
* ----------------------------------------------------------------------*/
|
||||
double ComputePressureUef::compute_scalar()
|
||||
{
|
||||
|
||||
temperature->compute_scalar();
|
||||
// if all pressures are external the scalar is found as normal
|
||||
if (ext_flags[0] && ext_flags[1] && ext_flags[2])
|
||||
return ComputePressure::compute_scalar();
|
||||
|
||||
// otherwise compute the full tensor and average desired components
|
||||
compute_vector();
|
||||
addstep(update->ntimestep+1);
|
||||
|
||||
int k =0;
|
||||
scalar = 0;
|
||||
if (ext_flags[0]) {
|
||||
scalar += vector[0];
|
||||
k++;
|
||||
}
|
||||
if (ext_flags[1]) {
|
||||
scalar += vector[1];
|
||||
k++;
|
||||
}
|
||||
if (ext_flags[2]) {
|
||||
scalar += vector[2];
|
||||
k++;
|
||||
}
|
||||
|
||||
scalar /= k;
|
||||
return scalar;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Compute the pressure tensor in the rotated coordinate system
|
||||
------------------------------------------------------------------------- */
|
||||
void ComputePressureUef::compute_vector()
|
||||
{
|
||||
invoked_vector = update->ntimestep;
|
||||
if (update->vflag_global != invoked_vector)
|
||||
error->all(FLERR,"Virial was not tallied on needed timestep");
|
||||
|
||||
if (force->kspace && kspace_virial && force->kspace->scalar_pressure_flag)
|
||||
error->all(FLERR,"Must use 'kspace_modify pressure/scalar no' for "
|
||||
"tensor components with kspace_style msm");
|
||||
|
||||
// invoke temperature if it hasn't been already
|
||||
|
||||
double *ke_tensor;
|
||||
if (keflag) {
|
||||
if (temperature->invoked_vector != update->ntimestep)
|
||||
temperature->compute_vector();
|
||||
ke_tensor = temperature->vector;
|
||||
}
|
||||
|
||||
if (dimension == 3) {
|
||||
inv_volume = 1.0 / (domain->xprd * domain->yprd * domain->zprd);
|
||||
virial_compute(6,3);
|
||||
if (in_fix)
|
||||
virial_rot(virial,rot);
|
||||
else
|
||||
{
|
||||
double r[3][3];
|
||||
( (FixNHUef*) modify->fix[ifix_uef])->get_rot(r);
|
||||
virial_rot(virial,r);
|
||||
}
|
||||
if (keflag) {
|
||||
for (int i = 0; i < 6; i++)
|
||||
vector[i] = (ke_tensor[i] + virial[i]) * inv_volume * nktv2p;
|
||||
} else
|
||||
for (int i = 0; i < 6; i++)
|
||||
vector[i] = virial[i] * inv_volume * nktv2p;
|
||||
} else {
|
||||
inv_volume = 1.0 / (domain->xprd * domain->yprd);
|
||||
virial_compute(4,2);
|
||||
if (keflag) {
|
||||
vector[0] = (ke_tensor[0] + virial[0]) * inv_volume * nktv2p;
|
||||
vector[1] = (ke_tensor[1] + virial[1]) * inv_volume * nktv2p;
|
||||
vector[3] = (ke_tensor[3] + virial[3]) * inv_volume * nktv2p;
|
||||
vector[2] = vector[4] = vector[5] = 0.0;
|
||||
} else {
|
||||
vector[0] = virial[0] * inv_volume * nktv2p;
|
||||
vector[1] = virial[1] * inv_volume * nktv2p;
|
||||
vector[3] = virial[3] * inv_volume * nktv2p;
|
||||
vector[2] = vector[4] = vector[5] = 0.0;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* get the current rotation matrix and store it
|
||||
------------------------------------------------------------------------- */
|
||||
void ComputePressureUef::update_rot()
|
||||
{
|
||||
( (FixNHUef*) modify->fix[ifix_uef])->get_rot(rot);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Transform the pressure tensor to the rotated coordinate system
|
||||
[P]rot = Q.[P].Q^t
|
||||
------------------------------------------------------------------------- */
|
||||
void ComputePressureUef::virial_rot(double *x, const double r[3][3])
|
||||
{
|
||||
|
||||
double t[3][3];
|
||||
|
||||
// [00 10 20 ] [ 0 3 4 ] [00 01 02 ]
|
||||
// [01 11 21 ] [ 3 1 5 ] [10 11 12 ]
|
||||
// [02 12 22 ] [ 4 5 2 ] [20 21 22 ]
|
||||
|
||||
for (int k = 0; k<3; ++k)
|
||||
{
|
||||
t[0][k] = x[0]*r[0][k] + x[3]*r[1][k] + x[4]*r[2][k];
|
||||
t[1][k] = x[3]*r[0][k] + x[1]*r[1][k] + x[5]*r[2][k];
|
||||
t[2][k] = x[4]*r[0][k] + x[5]*r[1][k] + x[2]*r[2][k];
|
||||
}
|
||||
x[0] = r[0][0]*t[0][0] + r[1][0]*t[1][0] + r[2][0]*t[2][0];
|
||||
x[3] = r[0][0]*t[0][1] + r[1][0]*t[1][1] + r[2][0]*t[2][1];
|
||||
x[4] = r[0][0]*t[0][2] + r[1][0]*t[1][2] + r[2][0]*t[2][2];
|
||||
x[1] = r[0][1]*t[0][1] + r[1][1]*t[1][1] + r[2][1]*t[2][1];
|
||||
x[5] = r[0][1]*t[0][2] + r[1][1]*t[1][2] + r[2][1]*t[2][2];
|
||||
x[2] = r[0][2]*t[0][2] + r[1][2]*t[1][2] + r[2][2]*t[2][2];
|
||||
}
|
||||
66
src/USER-UEF/compute_pressure_uef.h
Normal file
66
src/USER-UEF/compute_pressure_uef.h
Normal file
@ -0,0 +1,66 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef COMPUTE_CLASS
|
||||
|
||||
ComputeStyle(pressure/uef,ComputePressureUef)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_COMPUTE_PRESSURE_UEF_H
|
||||
#define LMP_COMPUTE_PRESSURE_UEF_H
|
||||
|
||||
#include "compute_pressure.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class ComputePressureUef : public ComputePressure {
|
||||
public:
|
||||
ComputePressureUef(class LAMMPS *, int, char **);
|
||||
virtual ~ComputePressureUef(){}
|
||||
virtual void init();
|
||||
virtual void compute_vector();
|
||||
virtual double compute_scalar();
|
||||
void update_rot();
|
||||
bool in_fix; //true if this compute is used in fix/nvt/npt
|
||||
|
||||
|
||||
protected:
|
||||
bool ext_flags[3]; // true if used in average output pressure
|
||||
void virial_rot(double*,const double[3][3]);
|
||||
int ifix_uef;
|
||||
double rot[3][3];
|
||||
};
|
||||
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
This class inherits most of the warnings from ComputePressure. The
|
||||
only additions are:
|
||||
|
||||
E: Can't use compute pressure/uef without defining a fix nvt/npt/uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
W: The temperature used in compute pressure/uef is not of style temp/uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
*/
|
||||
106
src/USER-UEF/compute_temp_uef.cpp
Normal file
106
src/USER-UEF/compute_temp_uef.cpp
Normal file
@ -0,0 +1,106 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include <string.h>
|
||||
#include <stdlib.h>
|
||||
#include "compute_temp_uef.h"
|
||||
#include "fix_nh_uef.h"
|
||||
#include "update.h"
|
||||
#include "modify.h"
|
||||
#include "fix.h"
|
||||
#include "error.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Base constructor initialized to use rotation matrix
|
||||
* ----------------------------------------------------------------------*/
|
||||
ComputeTempUef::ComputeTempUef(LAMMPS *lmp, int narg, char **arg) :
|
||||
ComputeTemp(lmp, narg, arg)
|
||||
{
|
||||
rot_flag=true;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Check for the uef fix
|
||||
* ----------------------------------------------------------------------*/
|
||||
void ComputeTempUef::init()
|
||||
{
|
||||
ComputeTemp::init();
|
||||
// check to make sure the other uef fix is on
|
||||
// borrowed from Pieter's nvt/sllod code
|
||||
int i=0;
|
||||
for (i=0; i<modify->nfix; i++) {
|
||||
if (strcmp(modify->fix[i]->style,"nvt/uef")==0)
|
||||
break;
|
||||
if (strcmp(modify->fix[i]->style,"npt/uef")==0)
|
||||
break;
|
||||
}
|
||||
if (i==modify->nfix)
|
||||
error->all(FLERR,"Can't use compute temp/uef without defining a fix nvt/npt/uef");
|
||||
ifix_uef=i;
|
||||
}
|
||||
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Compute the ke tensor in the proper coordinate system
|
||||
------------------------------------------------------------------------- */
|
||||
void ComputeTempUef::compute_vector()
|
||||
{
|
||||
ComputeTemp::compute_vector();
|
||||
if (rot_flag) {
|
||||
double rot[3][3];
|
||||
( (FixNHUef*) modify->fix[ifix_uef])->get_rot(rot);
|
||||
virial_rot(vector,rot);
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* turn the rotation matrix on or off to properly account for the
|
||||
* coordinate system of the velocities
|
||||
------------------------------------------------------------------------- */
|
||||
void ComputeTempUef::yes_rot()
|
||||
{
|
||||
rot_flag =true;
|
||||
}
|
||||
void ComputeTempUef::no_rot()
|
||||
{
|
||||
rot_flag =false;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Transform the pressure tensor to the rotated coordinate system
|
||||
[P]rot = Q.[P].Q^t
|
||||
------------------------------------------------------------------------- */
|
||||
void ComputeTempUef::virial_rot(double *x, const double r[3][3])
|
||||
{
|
||||
|
||||
double t[3][3];
|
||||
// [00 10 20 ] [ 0 3 4 ] [00 01 02 ]
|
||||
// [01 11 21 ] [ 3 1 5 ] [10 11 12 ]
|
||||
// [02 12 22 ] [ 4 5 2 ] [20 21 22 ]
|
||||
for (int k = 0; k<3; ++k) {
|
||||
t[0][k] = x[0]*r[0][k] + x[3]*r[1][k] + x[4]*r[2][k];
|
||||
t[1][k] = x[3]*r[0][k] + x[1]*r[1][k] + x[5]*r[2][k];
|
||||
t[2][k] = x[4]*r[0][k] + x[5]*r[1][k] + x[2]*r[2][k];
|
||||
}
|
||||
x[0] = r[0][0]*t[0][0] + r[1][0]*t[1][0] + r[2][0]*t[2][0];
|
||||
x[3] = r[0][0]*t[0][1] + r[1][0]*t[1][1] + r[2][0]*t[2][1];
|
||||
x[4] = r[0][0]*t[0][2] + r[1][0]*t[1][2] + r[2][0]*t[2][2];
|
||||
x[1] = r[0][1]*t[0][1] + r[1][1]*t[1][1] + r[2][1]*t[2][1];
|
||||
x[5] = r[0][1]*t[0][2] + r[1][1]*t[1][2] + r[2][1]*t[2][2];
|
||||
x[2] = r[0][2]*t[0][2] + r[1][2]*t[1][2] + r[2][2]*t[2][2];
|
||||
}
|
||||
60
src/USER-UEF/compute_temp_uef.h
Normal file
60
src/USER-UEF/compute_temp_uef.h
Normal file
@ -0,0 +1,60 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef COMPUTE_CLASS
|
||||
|
||||
ComputeStyle(temp/uef,ComputeTempUef)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_COMPUTE_TEMP_UEF_H
|
||||
#define LMP_COMPUTE_TEMP_UEF_H
|
||||
|
||||
#include "compute_temp.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class ComputeTempUef : public ComputeTemp {
|
||||
public:
|
||||
ComputeTempUef(class LAMMPS *, int, char **);
|
||||
virtual ~ComputeTempUef(){}
|
||||
virtual void init();
|
||||
virtual void compute_vector();
|
||||
void yes_rot();
|
||||
void no_rot();
|
||||
|
||||
|
||||
protected:
|
||||
bool rot_flag;
|
||||
void virial_rot(double*,const double[3][3]);
|
||||
int ifix_uef;
|
||||
};
|
||||
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
This class inherits most of the warnings from ComputePressure. The
|
||||
only addition is:
|
||||
|
||||
E: Can't use compute temp/uef without defining a fix nvt/npt/uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
*/
|
||||
114
src/USER-UEF/dump_cfg_uef.cpp
Normal file
114
src/USER-UEF/dump_cfg_uef.cpp
Normal file
@ -0,0 +1,114 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing Author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
|
||||
#include <math.h>
|
||||
#include <stdlib.h>
|
||||
#include <string.h>
|
||||
#include "dump_cfg.h"
|
||||
#include "atom.h"
|
||||
#include "domain.h"
|
||||
#include "modify.h"
|
||||
#include "compute.h"
|
||||
#include "fix.h"
|
||||
#include "error.h"
|
||||
#include "uef_utils.h"
|
||||
#include "dump_cfg_uef.h"
|
||||
#include "fix_nh_uef.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
enum{INT,DOUBLE,STRING,BIGINT}; // same as in DumpCustom
|
||||
|
||||
#define UNWRAPEXPAND 10.0
|
||||
#define ONEFIELD 32
|
||||
#define DELTA 1048576
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* base method is mostly fine, just need to find the FixNHUef
|
||||
* ----------------------------------------------------------------------*/
|
||||
void DumpCFGUef::init_style()
|
||||
{
|
||||
DumpCFG::init_style();
|
||||
|
||||
// check to make sure the other uef fix is on
|
||||
// borrowed from Pieter's nvt/sllod code
|
||||
int i=0;
|
||||
for (i=0; i<modify->nfix; i++)
|
||||
{
|
||||
if (strcmp(modify->fix[i]->style,"nvt/uef")==0)
|
||||
break;
|
||||
if (strcmp(modify->fix[i]->style,"npt/uef")==0)
|
||||
break;
|
||||
}
|
||||
if (i==modify->nfix)
|
||||
error->all(FLERR,"Can't use dump cfg/uef without defining a fix nvt/npt/uef");
|
||||
ifix_uef=i;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* this is really the only difference between the base class and this one.
|
||||
* since the output is in scaled coordinates, changing the simulation box
|
||||
* edges to the flow frame will put coordinates in the flow frame too.
|
||||
* ----------------------------------------------------------------------*/
|
||||
|
||||
void DumpCFGUef::write_header(bigint n)
|
||||
{
|
||||
// set scale factor used by AtomEye for CFG viz
|
||||
// default = 1.0
|
||||
// for peridynamics, set to pre-computed PD scale factor
|
||||
// so PD particles mimic C atoms
|
||||
// for unwrapped coords, set to UNWRAPEXPAND (10.0)
|
||||
// so molecules are not split across periodic box boundaries
|
||||
|
||||
double box[3][3],rot[3][3];
|
||||
((FixNHUef*) modify->fix[ifix_uef])->get_box(box);
|
||||
((FixNHUef*) modify->fix[ifix_uef])->get_rot(rot);
|
||||
// rot goes from "lab frame" to "upper triangular frame"
|
||||
// it's transpose takes the simulation box to the flow frame
|
||||
for (int i=0;i<3;i++)
|
||||
for(int j=i+1;j<3;j++)
|
||||
{
|
||||
double t=rot[i][j];
|
||||
rot[i][j]=rot[j][i];
|
||||
rot[j][i]=t;
|
||||
}
|
||||
UEF_utils::mul_m2(rot,box);
|
||||
|
||||
|
||||
double scale = 1.0;
|
||||
if (atom->peri_flag) scale = atom->pdscale;
|
||||
else if (unwrapflag == 1) scale = UNWRAPEXPAND;
|
||||
|
||||
char str[64];
|
||||
sprintf(str,"Number of particles = %s\n",BIGINT_FORMAT);
|
||||
fprintf(fp,str,n);
|
||||
fprintf(fp,"A = %g Angstrom (basic length-scale)\n",scale);
|
||||
// in box[][] columns are cell edges
|
||||
// in H0, rows are cell edges
|
||||
fprintf(fp,"H0(1,1) = %g A\n",box[0][0]);
|
||||
fprintf(fp,"H0(1,2) = %g A\n",box[1][0]);
|
||||
fprintf(fp,"H0(1,3) = %g A\n",box[2][0]);
|
||||
fprintf(fp,"H0(2,1) = %g A\n",box[0][1]);
|
||||
fprintf(fp,"H0(2,2) = %g A\n",box[1][1]);
|
||||
fprintf(fp,"H0(2,3) = %g A\n",box[2][1]);
|
||||
fprintf(fp,"H0(3,1) = %g A\n",box[0][2]);
|
||||
fprintf(fp,"H0(3,2) = %g A\n",box[1][2]);
|
||||
fprintf(fp,"H0(3,3) = %g A\n",box[2][2]);
|
||||
fprintf(fp,".NO_VELOCITY.\n");
|
||||
fprintf(fp,"entry_count = %d\n",nfield-2);
|
||||
for (int i = 0; i < nfield-5; i++)
|
||||
fprintf(fp,"auxiliary[%d] = %s\n",i,auxname[i]);
|
||||
}
|
||||
51
src/USER-UEF/dump_cfg_uef.h
Normal file
51
src/USER-UEF/dump_cfg_uef.h
Normal file
@ -0,0 +1,51 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing Author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef DUMP_CLASS
|
||||
|
||||
DumpStyle(cfg/uef,DumpCFGUef)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_DUMP_CFG_UEF_H
|
||||
#define LMP_DUMP_CFG_UEF_H
|
||||
|
||||
#include "dump_cfg.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class DumpCFGUef : public DumpCFG {
|
||||
public:
|
||||
DumpCFGUef(LAMMPS *lmp, int narg, char **arg) :
|
||||
DumpCFG(lmp, narg, arg){}
|
||||
void init_style();
|
||||
void write_header(bigint);
|
||||
|
||||
protected:
|
||||
int ifix_uef;
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
E: Can't use dump cfg/uef without defining a fix nvt/npt/uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
*/
|
||||
820
src/USER-UEF/fix_nh_uef.cpp
Normal file
820
src/USER-UEF/fix_nh_uef.cpp
Normal file
@ -0,0 +1,820 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
www.cs.sandia.gov/~sjplimp/lammps.html
|
||||
Steve Plimpton, sjplimp@sandia.gov, Sandia National Laboratories
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include <string.h>
|
||||
#include <stdlib.h>
|
||||
#include <math.h>
|
||||
#include "fix_nh_uef.h"
|
||||
#include "atom.h"
|
||||
#include "force.h"
|
||||
#include "group.h"
|
||||
#include "comm.h"
|
||||
#include "citeme.h"
|
||||
#include "irregular.h"
|
||||
#include "modify.h"
|
||||
#include "compute.h"
|
||||
#include "kspace.h"
|
||||
#include "update.h"
|
||||
#include "domain.h"
|
||||
#include "error.h"
|
||||
#include "output.h"
|
||||
#include "timer.h"
|
||||
#include "neighbor.h"
|
||||
#include "compute_pressure_uef.h"
|
||||
#include "compute_temp_uef.h"
|
||||
#include "uef_utils.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
using namespace FixConst;
|
||||
|
||||
enum{ISO,ANISO,TRICLINIC};
|
||||
|
||||
// citation info
|
||||
|
||||
static const char cite_user_uef_package[] =
|
||||
"USER-UEF package:\n\n"
|
||||
"@Article{NicholsonRutledge16,\n"
|
||||
"author = {David A. Nicholson and Gregory C. Rutledge},\n"
|
||||
"title = {Molecular simulation of flow-enhanced nucleation in n-eicosane melts under steady shear and uniaxial extension},\n"
|
||||
"journal = {The Journal of Chemical Physics},\n"
|
||||
"volume = {145},\n"
|
||||
"number = {24},\n"
|
||||
"pages = {244903},\n"
|
||||
"year = {2016}\n"
|
||||
"}\n\n";
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Parse fix specific keywords, do some error checking, and initalize
|
||||
* temp/pressure fixes
|
||||
---------------------------------------------------------------------- */
|
||||
FixNHUef::FixNHUef(LAMMPS *lmp, int narg, char **arg) :
|
||||
FixNH(lmp, narg, arg), uefbox(NULL)
|
||||
{
|
||||
if (lmp->citeme) lmp->citeme->add(cite_user_uef_package);
|
||||
|
||||
//initialization
|
||||
|
||||
erate[0] = erate[1] = 0;
|
||||
|
||||
// default values
|
||||
|
||||
strain[0]=strain[1]= 0;
|
||||
ext_flags[0]=ext_flags[1]=ext_flags[2] = true;
|
||||
|
||||
// need to initialize these
|
||||
|
||||
omega_dot[0]=omega_dot[1]=omega_dot[2]=0;
|
||||
|
||||
// parse fix nh/uef specific options
|
||||
|
||||
bool erate_flag = false;
|
||||
int iarg = 3;
|
||||
|
||||
while (iarg <narg) {
|
||||
if (strcmp(arg[iarg],"erate")==0) {
|
||||
if (iarg+3 > narg) error->all(FLERR,"Illegal fix nvt/npt/uef command");
|
||||
erate[0] = force->numeric(FLERR,arg[iarg+1]);
|
||||
erate[1] = force->numeric(FLERR,arg[iarg+2]);
|
||||
erate_flag = true;
|
||||
iarg += 3;
|
||||
} else if (strcmp(arg[iarg],"strain")==0) {
|
||||
if (iarg+3 > narg) error->all(FLERR,"Illegal fix nvt/npt/uef command");
|
||||
strain[0] = force->numeric(FLERR,arg[iarg+1]);
|
||||
strain[1] = force->numeric(FLERR,arg[iarg+2]);
|
||||
iarg += 3;
|
||||
} else if (strcmp(arg[iarg],"ext")==0) {
|
||||
if (iarg+2 > narg) error->all(FLERR,"Illegal fix nvt/npt/uef command");
|
||||
if (strcmp(arg[iarg+1],"x")==0)
|
||||
ext_flags[1] = ext_flags[2] = false;
|
||||
else if (strcmp(arg[iarg+1],"y")==0)
|
||||
ext_flags[0] = ext_flags[2] = false;
|
||||
else if (strcmp(arg[iarg+1],"z")==0)
|
||||
ext_flags[0] = ext_flags[1] = false;
|
||||
else if (strcmp(arg[iarg+1],"xy")==0)
|
||||
ext_flags[2] = false;
|
||||
else if (strcmp(arg[iarg+1],"xz")==0)
|
||||
ext_flags[1] = false;
|
||||
else if (strcmp(arg[iarg+1],"yz")==0)
|
||||
ext_flags[0] = false;
|
||||
else if (strcmp(arg[iarg+1],"xyz")!=0)
|
||||
error->all(FLERR,"Illegal fix nvt/npt/uef command");
|
||||
|
||||
iarg += 2;
|
||||
} else {
|
||||
|
||||
// skip to next argument; argument check for unknown keywords is done in FixNH
|
||||
|
||||
++iarg;
|
||||
}
|
||||
}
|
||||
|
||||
if (!erate_flag)
|
||||
error->all(FLERR,"Keyword erate must be set for fix npt/npt/uef command");
|
||||
|
||||
if (mtchain_default_flag) mtchain=1;
|
||||
|
||||
if (!domain->triclinic)
|
||||
error->all(FLERR,"Simulation box must be triclinic for fix/nvt/npt/uef");
|
||||
|
||||
// check for conditions that impose a deviatoric stress
|
||||
|
||||
if (pstyle == TRICLINIC)
|
||||
error->all(FLERR,"Only normal stresses can be controlled with fix/nvt/npt/uef");
|
||||
double erate_tmp[3];
|
||||
erate_tmp[0]=erate[0];
|
||||
erate_tmp[1]=erate[1];
|
||||
erate_tmp[2]=-erate[0]-erate[1];
|
||||
|
||||
if (pstyle == ANISO) {
|
||||
if (!(ext_flags[0] & ext_flags[1] & ext_flags[2]))
|
||||
error->all(FLERR,"The ext keyword may only be used with iso pressure control");
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
if (p_flag[k] && p_flag[j]) {
|
||||
double tol = 1e-6;
|
||||
if ( !nearly_equal(p_start[k],p_start[j],tol)
|
||||
|| !nearly_equal(p_stop[k],p_stop[j],tol))
|
||||
error->all(FLERR,"All controlled stresses must have the same "
|
||||
"value in fix/nvt/npt/uef");
|
||||
if ( !nearly_equal(erate_tmp[k],erate_tmp[j],tol)
|
||||
|| !nearly_equal(erate_tmp[k],erate_tmp[j],tol))
|
||||
error->all(FLERR,"Dimensions with controlled stresses must have"\
|
||||
" same strain rate in fix/nvt/npt/uef");
|
||||
}
|
||||
}
|
||||
|
||||
// conditions that produce a deviatoric stress have already been eliminated.
|
||||
|
||||
deviatoric_flag=0;
|
||||
|
||||
// need pre_exchange and irregular migration
|
||||
|
||||
pre_exchange_flag = 1;
|
||||
irregular = new Irregular(lmp);
|
||||
|
||||
// flag that I change the box here (in case of nvt)
|
||||
|
||||
box_change_shape = 1;
|
||||
|
||||
// initialize the UEFBox class which computes the box at each step
|
||||
|
||||
uefbox = new UEF_utils::UEFBox();
|
||||
uefbox->set_strain(strain[0],strain[1]);
|
||||
|
||||
// reset fixedpoint to the stagnation point. I don't allow fixedpoint
|
||||
// to be set by the user.
|
||||
|
||||
fixedpoint[0] = domain->boxlo[0];
|
||||
fixedpoint[1] = domain->boxlo[1];
|
||||
fixedpoint[2] = domain->boxlo[2];
|
||||
|
||||
// Create temp and pressure computes for nh/uef
|
||||
|
||||
int n = strlen(id) + 6;
|
||||
id_temp = new char[n];
|
||||
strcpy(id_temp,id);
|
||||
strcat(id_temp,"_temp");
|
||||
char **newarg = new char*[3];
|
||||
newarg[0] = id_temp;
|
||||
newarg[1] = (char *) "all";
|
||||
newarg[2] = (char *) "temp/uef";
|
||||
modify->add_compute(3,newarg);
|
||||
delete [] newarg;
|
||||
tcomputeflag = 1;
|
||||
|
||||
n = strlen(id) + 7;
|
||||
id_press = new char[n];
|
||||
strcpy(id_press,id);
|
||||
strcat(id_press,"_press");
|
||||
newarg = new char*[4];
|
||||
newarg[0] = id_press;
|
||||
newarg[1] = (char *) "all";
|
||||
newarg[2] = (char *) "pressure/uef";
|
||||
newarg[3] = id_temp;
|
||||
modify->add_compute(4,newarg);
|
||||
delete [] newarg;
|
||||
pcomputeflag = 1;
|
||||
|
||||
nevery = 1;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Erase the UEFBox object and get rid of the pressure compute if the nvt
|
||||
* version is being used. Everything else will be done in base destructor
|
||||
* ---------------------------------------------------------------------- */
|
||||
FixNHUef::~FixNHUef()
|
||||
{
|
||||
delete uefbox;
|
||||
if (pcomputeflag && !pstat_flag)
|
||||
{
|
||||
modify->delete_compute(id_press);
|
||||
delete [] id_press;
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Make the end_of_step() routine callable
|
||||
* ---------------------------------------------------------------------- */
|
||||
int FixNHUef::setmask()
|
||||
{
|
||||
int mask = FixNH::setmask();
|
||||
mask |= END_OF_STEP;
|
||||
return mask;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Run FixNH::init() and do more error checking. Set the pressure
|
||||
* pointer in the case that the nvt version is used
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::init()
|
||||
{
|
||||
FixNH::init();
|
||||
|
||||
|
||||
// find conflict with fix/deform or other box chaging fixes
|
||||
for (int i=0; i < modify->nfix; i++)
|
||||
{
|
||||
if (strcmp(modify->fix[i]->id,id) != 0)
|
||||
if (modify->fix[i]->box_change_shape != 0)
|
||||
error->all(FLERR,"Can't use another fix which changes box shape with fix/nvt/npt/uef");
|
||||
}
|
||||
|
||||
|
||||
// this will make the pressure compute for nvt
|
||||
if (!pstat_flag)
|
||||
if (pcomputeflag) {
|
||||
int icomp = modify->find_compute(id_press);
|
||||
if (icomp<0)
|
||||
error->all(FLERR,"Pressure ID for fix/nvt/uef doesn't exist");
|
||||
pressure = modify->compute[icomp];
|
||||
|
||||
if (strcmp(pressure->style,"pressure/uef") != 0)
|
||||
error->all(FLERR,"Using fix nvt/npt/uef without a compute pressure/uef");
|
||||
}
|
||||
|
||||
if (strcmp(temperature->style,"temp/uef") != 0)
|
||||
error->all(FLERR,"Using fix nvt/npt/uef without a compute temp/uef");
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Run FixNH::setup() make sure the box is OK and set the rotation matrix
|
||||
* for the first step
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::setup(int j)
|
||||
{
|
||||
double box[3][3];
|
||||
double vol = domain->xprd * domain->yprd * domain->zprd;
|
||||
uefbox->get_box(box,vol);
|
||||
double tol = 1e-4;
|
||||
// ensure the box is ok for uef
|
||||
bool isok = true;
|
||||
isok &= nearly_equal(domain->h[0],box[0][0],tol);
|
||||
isok &= nearly_equal(domain->h[1],box[1][1],tol);
|
||||
isok &= nearly_equal(domain->h[2],box[2][2],tol);
|
||||
isok &= nearly_equal(domain->xy,box[0][1],tol);
|
||||
isok &= nearly_equal(domain->yz,box[1][2],tol);
|
||||
isok &= nearly_equal(domain->xz,box[0][2],tol);
|
||||
if (!isok)
|
||||
error->all(FLERR,"Initial box is not close enough to the expected uef box");
|
||||
|
||||
uefbox->get_rot(rot);
|
||||
((ComputeTempUef*) temperature)->yes_rot();
|
||||
((ComputePressureUef*) pressure)->in_fix = true;
|
||||
((ComputePressureUef*) pressure)->update_rot();
|
||||
FixNH::setup(j);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* rotate -> initial integration step -> rotate back
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::initial_integrate(int vflag)
|
||||
{
|
||||
inv_rotate_x(rot);
|
||||
inv_rotate_v(rot);
|
||||
inv_rotate_f(rot);
|
||||
((ComputeTempUef*) temperature)->no_rot();
|
||||
FixNH::initial_integrate(vflag);
|
||||
rotate_x(rot);
|
||||
rotate_v(rot);
|
||||
rotate_f(rot);
|
||||
((ComputeTempUef*) temperature)->yes_rot();
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* rotate -> initial integration step -> rotate back (RESPA)
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::initial_integrate_respa(int vflag, int ilevel, int iloop)
|
||||
{
|
||||
inv_rotate_x(rot);
|
||||
inv_rotate_v(rot);
|
||||
inv_rotate_f(rot);
|
||||
((ComputeTempUef*) temperature)->no_rot();
|
||||
FixNH::initial_integrate_respa(vflag,ilevel,iloop);
|
||||
rotate_x(rot);
|
||||
rotate_v(rot);
|
||||
rotate_f(rot);
|
||||
((ComputeTempUef*) temperature)->yes_rot();
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* rotate -> final integration step -> rotate back
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::final_integrate()
|
||||
{
|
||||
// update rot here since it must directly follow the virial calculation
|
||||
((ComputePressureUef*) pressure)->update_rot();
|
||||
inv_rotate_v(rot);
|
||||
inv_rotate_f(rot);
|
||||
((ComputeTempUef*) temperature)->no_rot();
|
||||
FixNH::final_integrate();
|
||||
rotate_v(rot);
|
||||
rotate_f(rot);
|
||||
((ComputeTempUef*) temperature)->yes_rot();
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* at outer level: call this->final_integrate()
|
||||
* at other levels: rotate -> 2nd verlet step -> rotate back
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::final_integrate_respa(int ilevel, int iloop)
|
||||
{
|
||||
// set timesteps by level
|
||||
dtf = 0.5 * step_respa[ilevel] * force->ftm2v;
|
||||
dthalf = 0.5 * step_respa[ilevel];
|
||||
// outermost level - update eta_dot and omega_dot, apply via final_integrate
|
||||
// all other levels - NVE update of v
|
||||
if (ilevel == nlevels_respa-1) final_integrate();
|
||||
else
|
||||
{
|
||||
inv_rotate_v(rot);
|
||||
inv_rotate_f(rot);
|
||||
nve_v();
|
||||
rotate_v(rot);
|
||||
rotate_f(rot);
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
SLLOD velocity update in time-reversible (i think) increments
|
||||
v -> exp(-edot*dt/2)*v
|
||||
v -> v +f/m*dt
|
||||
v -> exp(-edot*dt/2)*v
|
||||
-----------------------------------------------------------------------*/
|
||||
void FixNHUef::nve_v()
|
||||
{
|
||||
double dtfm;
|
||||
double **v = atom->v;
|
||||
double **f = atom->f;
|
||||
double *rmass = atom->rmass;
|
||||
double *mass = atom->mass;
|
||||
int *type = atom->type;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
double ex = erate[0]*dtf/2;
|
||||
double ey = erate[1]*dtf/2;
|
||||
double ez = -ex-ey;
|
||||
double e0 = exp(-ex);
|
||||
double e1 = exp(-ey);
|
||||
double e2 = exp(-ez);
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
|
||||
if (rmass) {
|
||||
for (int i = 0; i < nlocal; i++) {
|
||||
if (mask[i] & groupbit) {
|
||||
dtfm = dtf / rmass[i];
|
||||
v[i][0] *= e0;
|
||||
v[i][1] *= e1;
|
||||
v[i][2] *= e2;
|
||||
v[i][0] += dtfm*f[i][0];
|
||||
v[i][1] += dtfm*f[i][1];
|
||||
v[i][2] += dtfm*f[i][2];
|
||||
v[i][0] *= e0;
|
||||
v[i][1] *= e1;
|
||||
v[i][2] *= e2;
|
||||
}
|
||||
}
|
||||
} else {
|
||||
for (int i = 0; i < nlocal; i++) {
|
||||
if (mask[i] & groupbit) {
|
||||
dtfm = dtf / mass[type[i]];
|
||||
v[i][0] *= e0;
|
||||
v[i][1] *= e1;
|
||||
v[i][2] *= e2;
|
||||
v[i][0] += dtfm*f[i][0];
|
||||
v[i][1] += dtfm*f[i][1];
|
||||
v[i][2] += dtfm*f[i][2];
|
||||
v[i][0] *= e0;
|
||||
v[i][1] *= e1;
|
||||
v[i][2] *= e2;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Don't actually move atoms in remap(), just change the box
|
||||
-----------------------------------------------------------------------*/
|
||||
void FixNHUef::remap()
|
||||
{
|
||||
double vol = domain->xprd * domain->yprd * domain->zprd;
|
||||
double domega = dto*(omega_dot[0]+omega_dot[1]+omega_dot[2])/3.;
|
||||
|
||||
// constant volume strain associated with barostat
|
||||
// box scaling
|
||||
double ex = dto*omega_dot[0]-domega;
|
||||
double ey = dto*omega_dot[1]-domega;
|
||||
uefbox->step_deform(ex,ey);
|
||||
strain[0] += ex;
|
||||
strain[1] += ey;
|
||||
|
||||
// volume change
|
||||
vol = vol*exp(3*domega);
|
||||
double box[3][3];
|
||||
uefbox->get_box(box,vol);
|
||||
domain->boxhi[0] = domain->boxlo[0]+box[0][0];
|
||||
domain->boxhi[1] = domain->boxlo[1]+box[1][1];
|
||||
domain->boxhi[2] = domain->boxlo[2]+box[2][2];
|
||||
domain->xy = box[0][1];
|
||||
domain->xz = box[0][2];
|
||||
domain->yz = box[1][2];
|
||||
domain->set_global_box();
|
||||
domain->set_local_box();
|
||||
uefbox->get_rot(rot);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
SLLOD position update in time-reversible (i think) increments
|
||||
x -> exp(edot*dt/2)*x
|
||||
x -> x + v*dt
|
||||
x -> exp(edot*dt/2)*x
|
||||
-----------------------------------------------------------------------*/
|
||||
void FixNHUef::nve_x()
|
||||
{
|
||||
double **x = atom->x;
|
||||
double **v = atom->v;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
double ex = erate[0]*dtv;
|
||||
strain[0] += ex;
|
||||
double e0 = exp((ex+omega_dot[0]*dtv)/2);
|
||||
double ey = erate[1]*dtv;
|
||||
strain[1] += ey;
|
||||
double e1 = exp((ey+omega_dot[1]*dtv)/2.);
|
||||
double ez = -ex -ey;
|
||||
double e2 = exp((ez+omega_dot[2]*dtv)/2.);
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
|
||||
// x update by full step only for atoms in group
|
||||
for (int i = 0; i < nlocal; i++) {
|
||||
if (mask[i] & groupbit) {
|
||||
x[i][0] *= e0;
|
||||
x[i][1] *= e1;
|
||||
x[i][2] *= e2;
|
||||
x[i][0] += dtv * v[i][0];
|
||||
x[i][1] += dtv * v[i][1];
|
||||
x[i][2] += dtv * v[i][2];
|
||||
x[i][0] *= e0;
|
||||
x[i][1] *= e1;
|
||||
x[i][2] *= e2;
|
||||
}
|
||||
}
|
||||
uefbox->step_deform(ex,ey);
|
||||
double box[3][3];
|
||||
double vol = domain->xprd * domain->yprd * domain->zprd;
|
||||
uefbox->get_box(box,vol);
|
||||
domain->boxhi[0] = domain->boxlo[0]+box[0][0];
|
||||
domain->boxhi[1] = domain->boxlo[1]+box[1][1];
|
||||
domain->boxhi[2] = domain->boxlo[2]+box[2][2];
|
||||
domain->xy = box[0][1];
|
||||
domain->xz = box[0][2];
|
||||
domain->yz = box[1][2];
|
||||
domain->set_global_box();
|
||||
domain->set_local_box();
|
||||
uefbox->get_rot(rot);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Do the lattice reduction if necessary.
|
||||
-----------------------------------------------------------------------*/
|
||||
void FixNHUef::pre_exchange()
|
||||
{
|
||||
// only need to reset things if the lattice needs to be reduced
|
||||
if (uefbox->reduce())
|
||||
{
|
||||
// go to lab frame
|
||||
inv_rotate_x(rot);
|
||||
inv_rotate_v(rot);
|
||||
inv_rotate_f(rot);
|
||||
// get & set the new box and rotation matrix
|
||||
double vol = domain->xprd * domain->yprd * domain->zprd;
|
||||
double box[3][3];
|
||||
uefbox->get_box(box,vol);
|
||||
domain->boxhi[0] = domain->boxlo[0]+box[0][0];
|
||||
domain->boxhi[1] = domain->boxlo[1]+box[1][1];
|
||||
domain->boxhi[2] = domain->boxlo[2]+box[2][2];
|
||||
domain->xy = box[0][1];
|
||||
domain->xz = box[0][2];
|
||||
domain->yz = box[1][2];
|
||||
domain->set_global_box();
|
||||
domain->set_local_box();
|
||||
uefbox->get_rot(rot);
|
||||
|
||||
// rotate to the new upper triangular frame
|
||||
rotate_v(rot);
|
||||
rotate_x(rot);
|
||||
rotate_f(rot);
|
||||
|
||||
// put all atoms in the new box
|
||||
double **x = atom->x;
|
||||
imageint *image = atom->image;
|
||||
int nlocal = atom->nlocal;
|
||||
for (int i=0; i<nlocal; i++) domain->remap(x[i],image[i]);
|
||||
|
||||
// move atoms to the right processors
|
||||
domain->x2lamda(atom->nlocal);
|
||||
irregular->migrate_atoms();
|
||||
domain->lamda2x(atom->nlocal);
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* The following are routines to rotate between the lab and upper triangular
|
||||
* (UT) frames. For most of the time the simulation is in the UT frame.
|
||||
* To get to the lab frame, apply the inv_rotate_[..](rot) and to
|
||||
* get back to the UT frame apply rotate_[..](rot).
|
||||
*
|
||||
* Note: the rotate_x() functions also apply a shift to/from the fixedpoint
|
||||
* to make the integration a little simpler.
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::rotate_x(double r[3][3])
|
||||
{
|
||||
double **x = atom->x;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
|
||||
double xn[3];
|
||||
for (int i=0;i<nlocal;i++)
|
||||
{
|
||||
if (mask[i] & groupbit)
|
||||
{
|
||||
xn[0]=r[0][0]*x[i][0]+r[0][1]*x[i][1]+r[0][2]*x[i][2];
|
||||
xn[1]=r[1][0]*x[i][0]+r[1][1]*x[i][1]+r[1][2]*x[i][2];
|
||||
xn[2]=r[2][0]*x[i][0]+r[2][1]*x[i][1]+r[2][2]*x[i][2];
|
||||
x[i][0]=xn[0]+domain->boxlo[0];
|
||||
x[i][1]=xn[1]+domain->boxlo[1];
|
||||
x[i][2]=xn[2]+domain->boxlo[2];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void FixNHUef::inv_rotate_x(double r[3][3])
|
||||
{
|
||||
double **x = atom->x;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
|
||||
double xn[3];
|
||||
for (int i=0;i<nlocal;i++)
|
||||
{
|
||||
if (mask[i] & groupbit)
|
||||
{
|
||||
x[i][0] -= domain->boxlo[0];
|
||||
x[i][1] -= domain->boxlo[1];
|
||||
x[i][2] -= domain->boxlo[2];
|
||||
xn[0]=r[0][0]*x[i][0]+r[1][0]*x[i][1]+r[2][0]*x[i][2];
|
||||
xn[1]=r[0][1]*x[i][0]+r[1][1]*x[i][1]+r[2][1]*x[i][2];
|
||||
xn[2]=r[0][2]*x[i][0]+r[1][2]*x[i][1]+r[2][2]*x[i][2];
|
||||
x[i][0]=xn[0];
|
||||
x[i][1]=xn[1];
|
||||
x[i][2]=xn[2];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void FixNHUef::rotate_v(double r[3][3])
|
||||
{
|
||||
double **v = atom->v;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
|
||||
double vn[3];
|
||||
for (int i=0;i<nlocal;i++)
|
||||
{
|
||||
if (mask[i] & groupbit)
|
||||
{
|
||||
vn[0]=r[0][0]*v[i][0]+r[0][1]*v[i][1]+r[0][2]*v[i][2];
|
||||
vn[1]=r[1][0]*v[i][0]+r[1][1]*v[i][1]+r[1][2]*v[i][2];
|
||||
vn[2]=r[2][0]*v[i][0]+r[2][1]*v[i][1]+r[2][2]*v[i][2];
|
||||
v[i][0]=vn[0]; v[i][1]=vn[1]; v[i][2]=vn[2];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void FixNHUef::inv_rotate_v(double r[3][3])
|
||||
{
|
||||
double **v = atom->v;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
|
||||
double vn[3];
|
||||
for (int i=0;i<nlocal;i++)
|
||||
{
|
||||
if (mask[i] & groupbit)
|
||||
{
|
||||
vn[0]=r[0][0]*v[i][0]+r[1][0]*v[i][1]+r[2][0]*v[i][2];
|
||||
vn[1]=r[0][1]*v[i][0]+r[1][1]*v[i][1]+r[2][1]*v[i][2];
|
||||
vn[2]=r[0][2]*v[i][0]+r[1][2]*v[i][1]+r[2][2]*v[i][2];
|
||||
v[i][0]=vn[0]; v[i][1]=vn[1]; v[i][2]=vn[2];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void FixNHUef::rotate_f(double r[3][3])
|
||||
{
|
||||
double **f = atom->f;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
|
||||
double fn[3];
|
||||
for (int i=0;i<nlocal;i++)
|
||||
{
|
||||
if (mask[i] & groupbit)
|
||||
{
|
||||
fn[0]=r[0][0]*f[i][0]+r[0][1]*f[i][1]+r[0][2]*f[i][2];
|
||||
fn[1]=r[1][0]*f[i][0]+r[1][1]*f[i][1]+r[1][2]*f[i][2];
|
||||
fn[2]=r[2][0]*f[i][0]+r[2][1]*f[i][1]+r[2][2]*f[i][2];
|
||||
f[i][0]=fn[0]; f[i][1]=fn[1]; f[i][2]=fn[2];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void FixNHUef::inv_rotate_f(double r[3][3])
|
||||
{
|
||||
double **f = atom->f;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
double fn[3];
|
||||
for (int i=0;i<nlocal;i++)
|
||||
{
|
||||
if (mask[i] & groupbit)
|
||||
{
|
||||
fn[0]=r[0][0]*f[i][0]+r[1][0]*f[i][1]+r[2][0]*f[i][2];
|
||||
fn[1]=r[0][1]*f[i][0]+r[1][1]*f[i][1]+r[2][1]*f[i][2];
|
||||
fn[2]=r[0][2]*f[i][0]+r[1][2]*f[i][1]+r[2][2]*f[i][2];
|
||||
f[i][0]=fn[0]; f[i][1]=fn[1]; f[i][2]=fn[2];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Increase the size of the restart list to add in the strains
|
||||
* ---------------------------------------------------------------------- */
|
||||
int FixNHUef::size_restart_global()
|
||||
{
|
||||
return FixNH::size_restart_global() +2;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* Pack the strains after packing the default FixNH values
|
||||
* ---------------------------------------------------------------------- */
|
||||
int FixNHUef::pack_restart_data(double *list)
|
||||
{
|
||||
int n = FixNH::pack_restart_data(list);
|
||||
list[n++] = strain[0];
|
||||
list[n++] = strain[1];
|
||||
return n;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* read and set the strains after the default FixNH values
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::restart(char *buf)
|
||||
{
|
||||
int n = size_restart_global();
|
||||
FixNH::restart(buf);
|
||||
double *list = (double *) buf;
|
||||
strain[0] = list[n-2];
|
||||
strain[1] = list[n-1];
|
||||
uefbox->set_strain(strain[0],strain[1]);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* If the step writes a restart, reduce the box beforehand. This makes sure
|
||||
* the unique box shape can be found once the restart is read and that
|
||||
* all of the atoms lie within the box.
|
||||
* This may only be necessary for RESPA runs, but I'm leaving it in anyway.
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::end_of_step()
|
||||
{
|
||||
if (update->ntimestep==output->next_restart)
|
||||
{
|
||||
pre_exchange();
|
||||
domain->x2lamda(atom->nlocal);
|
||||
domain->pbc();
|
||||
timer->stamp();
|
||||
comm->exchange();
|
||||
comm->borders();
|
||||
domain->lamda2x(atom->nlocal+atom->nghost);
|
||||
timer->stamp(Timer::COMM);
|
||||
neighbor->build();
|
||||
timer->stamp(Timer::NEIGH);
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* reduce the simulation box after a run is complete. otherwise it won't
|
||||
* be possible to resume from a write_restart since the initialization of
|
||||
* the simulation box requires reduced simulation box
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::post_run()
|
||||
{
|
||||
pre_exchange();
|
||||
domain->x2lamda(atom->nlocal);
|
||||
domain->pbc();
|
||||
timer->stamp();
|
||||
comm->exchange();
|
||||
comm->borders();
|
||||
domain->lamda2x(atom->nlocal+atom->nghost);
|
||||
timer->stamp(Timer::COMM);
|
||||
neighbor->build();
|
||||
timer->stamp(Timer::NEIGH);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* public read for rotation matrix
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::get_rot(double r[3][3])
|
||||
{
|
||||
r[0][0] = rot[0][0];
|
||||
r[0][1] = rot[0][1];
|
||||
r[0][2] = rot[0][2];
|
||||
r[1][0] = rot[1][0];
|
||||
r[1][1] = rot[1][1];
|
||||
r[1][2] = rot[1][2];
|
||||
r[2][0] = rot[2][0];
|
||||
r[2][1] = rot[2][1];
|
||||
r[2][2] = rot[2][2];
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* public read for ext flags
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::get_ext_flags(bool* e)
|
||||
{
|
||||
e[0] = ext_flags[0];
|
||||
e[1] = ext_flags[1];
|
||||
e[2] = ext_flags[2];
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* public read for simulation box
|
||||
* ---------------------------------------------------------------------- */
|
||||
void FixNHUef::get_box(double b[3][3])
|
||||
{
|
||||
double box[3][3];
|
||||
double vol = domain->xprd * domain->yprd * domain->zprd;
|
||||
uefbox->get_box(box,vol);
|
||||
b[0][0] = box[0][0];
|
||||
b[0][1] = box[0][1];
|
||||
b[0][2] = box[0][2];
|
||||
b[1][0] = box[1][0];
|
||||
b[1][1] = box[1][1];
|
||||
b[1][2] = box[1][2];
|
||||
b[2][0] = box[2][0];
|
||||
b[2][1] = box[2][1];
|
||||
b[2][2] = box[2][2];
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
* comparing floats
|
||||
* it's imperfect, but should work provided no infinities
|
||||
* ---------------------------------------------------------------------- */
|
||||
bool FixNHUef::nearly_equal(double a, double b, double epsilon)
|
||||
{
|
||||
double absa = fabs(a);
|
||||
double absb = fabs(b);
|
||||
double diff = fabs(a-b);
|
||||
if (a == b) return true;
|
||||
else if ( (absa+absb) < epsilon)
|
||||
return diff < epsilon*epsilon;
|
||||
else
|
||||
return diff/(absa+absb) < epsilon;
|
||||
}
|
||||
128
src/USER-UEF/fix_nh_uef.h
Normal file
128
src/USER-UEF/fix_nh_uef.h
Normal file
@ -0,0 +1,128 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
www.cs.sandia.gov/~sjplimp/lammps.html
|
||||
Steve Plimpton, sjplimp@sandia.gov, Sandia National Laboratories
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing Author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
|
||||
#ifndef LMP_FIX_NH_UEF_H
|
||||
#define LMP_FIX_NH_UEF_H
|
||||
|
||||
#include "fix_nh.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
// forward declaration
|
||||
namespace UEF_utils {
|
||||
class UEFBox;
|
||||
};
|
||||
|
||||
class FixNHUef : public FixNH {
|
||||
public:
|
||||
FixNHUef(class LAMMPS *, int, char **);
|
||||
virtual ~FixNHUef();
|
||||
virtual int setmask();
|
||||
virtual void init();
|
||||
virtual void setup(int);
|
||||
virtual void pre_exchange();
|
||||
virtual int pack_restart_data(double*);
|
||||
virtual void restart(char *);
|
||||
virtual void end_of_step();
|
||||
virtual void initial_integrate(int);
|
||||
virtual void final_integrate();
|
||||
virtual void initial_integrate_respa(int, int, int);
|
||||
virtual void final_integrate_respa(int, int);
|
||||
virtual void post_run();
|
||||
void get_rot(double[3][3]);
|
||||
void get_ext_flags(bool*);
|
||||
void get_box(double[3][3]);
|
||||
|
||||
protected:
|
||||
virtual void remap();
|
||||
virtual int size_restart_global();
|
||||
virtual void nve_x();
|
||||
virtual void nve_v();
|
||||
void rotate_x(double [3][3]);
|
||||
void inv_rotate_x(double[3][3]);
|
||||
void rotate_v(double[3][3]);
|
||||
void inv_rotate_v(double[3][3]);
|
||||
void rotate_f(double[3][3]);
|
||||
void inv_rotate_f(double[3][3]);
|
||||
double strain[2],erate[2]; // strain/strain rate : [e_x, e_y]
|
||||
// always assume traceless e_z = -e_x-e_y
|
||||
|
||||
int rem; //this is for the narg kluge
|
||||
|
||||
UEF_utils::UEFBox *uefbox; // interface for the special simulation box
|
||||
|
||||
double rot[3][3]; // rotation matrix
|
||||
bool ext_flags[3]; // flags for external "free surfaces"
|
||||
bool nearly_equal(double,double,double);
|
||||
//bool rotate_output; // experimental feature. Too many issues for now
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
This is a base class for FixNH so it will inherit most of its error/warning messages along with the following:
|
||||
|
||||
E: Illegal fix nvt/npt/uef command
|
||||
|
||||
Self-explanatory
|
||||
|
||||
E: Keyword erate must be set for fix nvt/npt/uef command
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Simulation box must be triclinic for fix/nvt/npt/uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Only normal stresses can be controlled with fix/nvt/npt/uef
|
||||
|
||||
The keywords xy xz and yz cannot be used for pressure control
|
||||
|
||||
E: The ext keyword may only be used with iso pressure control
|
||||
|
||||
Self-explanatory
|
||||
|
||||
E: All controlled stresses must have the same value in fix/nvt/npt/uef
|
||||
|
||||
Stress control is only possible when the stress specified for each dimension is the same
|
||||
|
||||
E: Dimensions with controlled stresses must have same strain rate in fix/nvt/npt/uef
|
||||
|
||||
Stress-controlled dimensions with the same strain rate must have the same target stress.
|
||||
|
||||
E: Can't use another fix which changes box shape with fix/nvt/npt/uef
|
||||
|
||||
The fix npt/nvt/uef command must have full control over the box shape. You cannot use a simultaneous fix deform command, for example.
|
||||
|
||||
E: Pressure ID for fix/nvt/uef doesn't exist
|
||||
|
||||
The compute pressure introduced via fix_modify does not exist
|
||||
|
||||
E: Using fix nvt/npt/uef without a compute pressure/uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Using fix nvt/npt/uef without a compute temp/uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Initial box is not close enough to the expected uef box
|
||||
|
||||
The initial box does not correspond to the shape required by the value of the strain keyword. If the default strain value of zero was used, the initial box is not cubic.
|
||||
|
||||
*/
|
||||
28
src/USER-UEF/fix_npt_uef.cpp
Normal file
28
src/USER-UEF/fix_npt_uef.cpp
Normal file
@ -0,0 +1,28 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include "fix_npt_uef.h"
|
||||
#include "error.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
FixNPTUef::FixNPTUef(LAMMPS *lmp, int narg, char **arg) :
|
||||
FixNHUef(lmp, narg, arg)
|
||||
{
|
||||
if (!tstat_flag)
|
||||
error->all(FLERR,"Temperature control must be used with fix npt/uef");
|
||||
if (!pstat_flag)
|
||||
error->all(FLERR,"Pressure control must be used with fix npt/uef");
|
||||
}
|
||||
50
src/USER-UEF/fix_npt_uef.h
Normal file
50
src/USER-UEF/fix_npt_uef.h
Normal file
@ -0,0 +1,50 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef FIX_CLASS
|
||||
|
||||
FixStyle(npt/uef,FixNPTUef)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_FIX_NPT_UEF_H
|
||||
#define LMP_FIX_NPT_UEF_H
|
||||
|
||||
#include "fix_nh_uef.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class FixNPTUef : public FixNHUef {
|
||||
public:
|
||||
FixNPTUef(class LAMMPS *, int, char **);
|
||||
~FixNPTUef() {}
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
E: Temperature control must be used with fix npt uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Pressure control must be used with fix npt uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
*/
|
||||
30
src/USER-UEF/fix_nvt_uef.cpp
Normal file
30
src/USER-UEF/fix_nvt_uef.cpp
Normal file
@ -0,0 +1,30 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include "fix_nvt_uef.h"
|
||||
#include "error.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
FixNVTUef::FixNVTUef(LAMMPS *lmp, int narg, char **arg) :
|
||||
FixNHUef(lmp, narg, arg)
|
||||
{
|
||||
if (!tstat_flag)
|
||||
error->all(FLERR,"Temperature control must be used with fix nvt/uef");
|
||||
if (pstat_flag)
|
||||
error->all(FLERR,"Pressure control can't be used with fix nvt/uef");
|
||||
}
|
||||
|
||||
|
||||
50
src/USER-UEF/fix_nvt_uef.h
Normal file
50
src/USER-UEF/fix_nvt_uef.h
Normal file
@ -0,0 +1,50 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef FIX_CLASS
|
||||
|
||||
FixStyle(nvt/uef,FixNVTUef)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_FIX_NVT_UEF_H
|
||||
#define LMP_FIX_NVT_UEF_H
|
||||
|
||||
#include "fix_nh_uef.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class FixNVTUef : public FixNHUef {
|
||||
public:
|
||||
FixNVTUef(class LAMMPS *, int, char **);
|
||||
~FixNVTUef(){}
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
E: Temperature control must be used with fix npt uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Pressure control must be used with fix npt uef
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
*/
|
||||
366
src/USER-UEF/uef_utils.cpp
Normal file
366
src/USER-UEF/uef_utils.cpp
Normal file
@ -0,0 +1,366 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
-------------------------------------------------------------------------
|
||||
|
||||
This class contains functions to calculate the evolution of the periodic
|
||||
simulation box under elongational flow as described by Matthew Dobson
|
||||
in the arXiv preprint at http://arxiv.org/abs/1408.7078
|
||||
|
||||
Additionally, there are methods to do a lattice reduction to further
|
||||
reduce the simulation box using the method of Igor Semaev at
|
||||
http://link.springer.com/chapter/10.1007%2F3-540-44670-2_13
|
||||
*/
|
||||
|
||||
#include <math.h>
|
||||
#include "uef_utils.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
namespace UEF_utils{
|
||||
|
||||
UEFBox::UEFBox()
|
||||
{
|
||||
// initial box (also an inverse eigenvector matrix of automorphisms)
|
||||
double x = 0.327985277605681;
|
||||
double y = 0.591009048506103;
|
||||
double z = 0.736976229099578;
|
||||
l0[0][0]= z; l0[0][1]= y; l0[0][2]= x;
|
||||
l0[1][0]=-x; l0[1][1]= z; l0[1][2]=-y;
|
||||
l0[2][0]=-y; l0[2][1]= x; l0[2][2]= z;
|
||||
// spectra of the two automorpisms (log of eigenvalues)
|
||||
w1[0]=-1.177725211523360;
|
||||
w1[1]=-0.441448620566067;
|
||||
w1[2]= 1.619173832089425;
|
||||
w2[0]= w1[1];
|
||||
w2[1]= w1[2];
|
||||
w2[2]= w1[0];
|
||||
// initialize theta
|
||||
// strain = w1 * theta1 + w2 * theta2
|
||||
theta[0]=theta[1]=0;
|
||||
|
||||
|
||||
//set up the initial box l and change of basis matrix r
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
{
|
||||
l[k][j] = l0[k][j];
|
||||
r[j][k]=(j==k);
|
||||
}
|
||||
|
||||
// get the initial rotation and upper triangular matrix
|
||||
rotation_matrix(rot, lrot ,l);
|
||||
|
||||
// this is just a way to calculate the automorphisms
|
||||
// themselves, which play a minor role in the calculations
|
||||
// it's overkill, but only called once
|
||||
double t1[3][3];
|
||||
double t1i[3][3];
|
||||
double t2[3][3];
|
||||
double t2i[3][3];
|
||||
double l0t[3][3];
|
||||
for (int k=0; k<3; ++k)
|
||||
for (int j=0; j<3; ++j)
|
||||
{
|
||||
t1[k][j] = exp(w1[k])*l0[k][j];
|
||||
t1i[k][j] = exp(-w1[k])*l0[k][j];
|
||||
t2[k][j] = exp(w2[k])*l0[k][j];
|
||||
t2i[k][j] = exp(-w2[k])*l0[k][j];
|
||||
l0t[k][j] = l0[j][k];
|
||||
}
|
||||
mul_m2(l0t,t1);
|
||||
mul_m2(l0t,t1i);
|
||||
mul_m2(l0t,t2);
|
||||
mul_m2(l0t,t2i);
|
||||
for (int k=0; k<3; ++k)
|
||||
for (int j=0; j<3; ++j)
|
||||
{
|
||||
a1[k][j] = round(t1[k][j]);
|
||||
a1i[k][j] = round(t1i[k][j]);
|
||||
a2[k][j] = round(t2[k][j]);
|
||||
a2i[k][j] = round(t2i[k][j]);
|
||||
}
|
||||
|
||||
// winv used to transform between
|
||||
// strain increments and theta increments
|
||||
winv[0][0] = w2[1];
|
||||
winv[0][1] = -w2[0];
|
||||
winv[1][0] = -w1[1];
|
||||
winv[1][1] = w1[0];
|
||||
double d = w1[0]*w2[1] - w2[0]*w1[1];
|
||||
for (int k=0;k<2;k++)
|
||||
for (int j=0;j<2;j++)
|
||||
winv[k][j] /= d;
|
||||
}
|
||||
|
||||
// get volume-correct r basis in: basis*cbrt(vol) = q*r
|
||||
void UEFBox::get_box(double x[3][3], double v)
|
||||
{
|
||||
v = cbrtf(v);
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
x[k][j] = lrot[k][j]*v;
|
||||
}
|
||||
|
||||
// get rotation matrix q in: basis = q*r
|
||||
void UEFBox::get_rot(double x[3][3])
|
||||
{
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
x[k][j]=rot[k][j];
|
||||
}
|
||||
|
||||
// diagonal, incompressible deformation
|
||||
void UEFBox::step_deform(const double ex, const double ey)
|
||||
{
|
||||
// increment theta values used in the reduction
|
||||
theta[0] +=winv[0][0]*ex + winv[0][1]*ey;
|
||||
theta[1] +=winv[1][0]*ex + winv[1][1]*ey;
|
||||
|
||||
// deformation of the box. reduce() needs to
|
||||
// be called regularly or calculation will become
|
||||
// unstable
|
||||
double eps[3];
|
||||
eps[0]=ex; eps[1] = ey; eps[2] = -ex-ey;
|
||||
for (int k=0;k<3;k++)
|
||||
{
|
||||
eps[k] = exp(eps[k]);
|
||||
l[k][0] = eps[k]*l[k][0];
|
||||
l[k][1] = eps[k]*l[k][1];
|
||||
l[k][2] = eps[k]*l[k][2];
|
||||
}
|
||||
rotation_matrix(rot,lrot, l);
|
||||
}
|
||||
// reuduce the current basis
|
||||
bool UEFBox::reduce()
|
||||
{
|
||||
// determine how many times to apply the automorphisms
|
||||
// and find new theta values
|
||||
int f1 = round(theta[0]);
|
||||
int f2 = round(theta[1]);
|
||||
theta[0] -= f1;
|
||||
theta[1] -= f2;
|
||||
|
||||
// store old change or basis matrix to determine if it
|
||||
// changes
|
||||
int r0[3][3];
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
r0[k][j]=r[k][j];
|
||||
|
||||
// this modifies the old change basis matrix to
|
||||
// handle the case where the automorphism transforms
|
||||
// the box but the reduced basis doesn't change
|
||||
// (r0 should still equal r at the end)
|
||||
if (f1 > 0) for (int k=0;k<f1;k++) mul_m2 (a1,r0);
|
||||
if (f1 < 0) for (int k=0;k<-f1;k++) mul_m2 (a1i,r0);
|
||||
if (f2 > 0) for (int k=0;k<f2;k++) mul_m2 (a2,r0);
|
||||
if (f2 < 0) for (int k=0;k<-f2;k++) mul_m2 (a2i,r0);
|
||||
|
||||
// robust reduction to the box defined by Dobson
|
||||
for (int k=0;k<3;k++)
|
||||
{
|
||||
double eps = exp(theta[0]*w1[k]+theta[1]*w2[k]);
|
||||
l[k][0] = eps*l0[k][0];
|
||||
l[k][1] = eps*l0[k][1];
|
||||
l[k][2] = eps*l0[k][2];
|
||||
}
|
||||
// further reduce the box using greedy reduction and check
|
||||
// if it changed from the last step using the change of basis
|
||||
// matrices r and r0
|
||||
greedy(l,r);
|
||||
rotation_matrix(rot,lrot, l);
|
||||
return !mat_same(r,r0);
|
||||
}
|
||||
void UEFBox::set_strain(const double ex, const double ey)
|
||||
{
|
||||
theta[0] =winv[0][0]*ex + winv[0][1]*ey;
|
||||
theta[1] =winv[1][0]*ex + winv[1][1]*ey;
|
||||
theta[0] -= round(theta[0]);
|
||||
theta[1] -= round(theta[1]);
|
||||
|
||||
for (int k=0;k<3;k++)
|
||||
{
|
||||
double eps = exp(theta[0]*w1[k]+theta[1]*w2[k]);
|
||||
l[k][0] = eps*l0[k][0];
|
||||
l[k][1] = eps*l0[k][1];
|
||||
l[k][2] = eps*l0[k][2];
|
||||
}
|
||||
greedy(l,r);
|
||||
rotation_matrix(rot,lrot, l);
|
||||
}
|
||||
|
||||
// this is just qr reduction using householder reflections
|
||||
// m is input matrix, q is a rotation, r is upper triangular
|
||||
// q*m = r
|
||||
void rotation_matrix(double q[3][3], double r[3][3], const double m[3][3])
|
||||
{
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
r[k][j] = m[k][j];
|
||||
|
||||
double a = -sqrt(col_prod(r,0,0))*r[0][0]/fabs(r[0][0]);
|
||||
double v[3];
|
||||
v[0] = r[0][0]-a;
|
||||
v[1] = r[1][0];
|
||||
v[2] = r[2][0];
|
||||
a = sqrt(v[0]*v[0]+v[1]*v[1]+v[2]*v[2]);
|
||||
v[0] /= a; v[1] /= a; v[2] /= a;
|
||||
double qt[3][3];
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
{
|
||||
qt[k][j] = (k==j) - 2*v[k]*v[j];
|
||||
q[k][j]= qt[k][j];
|
||||
}
|
||||
mul_m2(qt,r);
|
||||
a = -sqrt(r[1][1]*r[1][1] + r[2][1]*r[2][1])*r[1][1]/fabs(r[1][1]);
|
||||
v[0] = 0;
|
||||
v[1] = r[1][1] - a;
|
||||
v[2] = r[2][1];
|
||||
a = sqrt(v[1]*v[1]+v[2]*v[2]);
|
||||
v[1] /= a;
|
||||
v[2] /= a;
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
qt[k][j] = (k==j) - 2*v[k]*v[j];
|
||||
mul_m2(qt,r);
|
||||
mul_m2(qt,q);
|
||||
// this makes r have positive diagonals
|
||||
// q*m = r <==> (-q)*m = (-r) will hold row-wise
|
||||
if (r[0][0] < 0){ neg_row(q,0); neg_row(r,0); }
|
||||
if (r[1][1] < 0){ neg_row(q,1); neg_row(r,1); }
|
||||
if (r[2][2] < 0){ neg_row(q,2); neg_row(r,2); }
|
||||
}
|
||||
|
||||
|
||||
|
||||
//sort columns in order of increasing length
|
||||
void col_sort(double b[3][3],int r[3][3])
|
||||
{
|
||||
if (col_prod(b,0,0)>col_prod(b,1,1))
|
||||
{
|
||||
col_swap(b,0,1);
|
||||
col_swap(r,0,1);
|
||||
}
|
||||
if (col_prod(b,0,0)>col_prod(b,2,2))
|
||||
{
|
||||
col_swap(b,0,2);
|
||||
col_swap(r,0,2);
|
||||
}
|
||||
if (col_prod(b,1,1)>col_prod(b,2,2))
|
||||
{
|
||||
col_swap(b,1,2);
|
||||
col_swap(r,1,2);
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
// 1-2 reduction (Graham-Schmidt)
|
||||
void red12(double b[3][3],int r[3][3])
|
||||
{
|
||||
int y = round(col_prod(b,0,1)/col_prod(b,0,0));
|
||||
b[0][1] -= y*b[0][0];
|
||||
b[1][1] -= y*b[1][0];
|
||||
b[2][1] -= y*b[2][0];
|
||||
|
||||
r[0][1] -= y*r[0][0];
|
||||
r[1][1] -= y*r[1][0];
|
||||
r[2][1] -= y*r[2][0];
|
||||
if (col_prod(b,1,1) < col_prod(b,0,0))
|
||||
{
|
||||
col_swap(b,0,1);
|
||||
col_swap(r,0,1);
|
||||
red12(b,r);
|
||||
}
|
||||
}
|
||||
|
||||
// The Semaev condition for a 3-reduced basis
|
||||
void red3(double b[3][3], int r[3][3])
|
||||
{
|
||||
double b11 = col_prod(b,0,0);
|
||||
double b22 = col_prod(b,1,1);
|
||||
double b12 = col_prod(b,0,1);
|
||||
double b13 = col_prod(b,0,2);
|
||||
double b23 = col_prod(b,1,2);
|
||||
|
||||
double y2 =-(b23/b22-b12/b22*b13/b11)/(1-b12/b11*b12/b22);
|
||||
double y1 =-(b13/b11-b12/b11*b23/b22)/(1-b12/b11*b12/b22);
|
||||
|
||||
int x1=0;
|
||||
int x2=0;
|
||||
double min = col_prod(b,2,2);
|
||||
int x1v[2];
|
||||
int x2v[2];
|
||||
x1v[0] = floor(y1); x1v[1] = x1v[0]+1;
|
||||
x2v[0] = floor(y2); x2v[1] = x2v[0]+1;
|
||||
for (int k=0;k<2;k++)
|
||||
for (int j=0;j<2;j++)
|
||||
{
|
||||
double a[3];
|
||||
a[0] = b[0][2] + x1v[k]*b[0][0] + x2v[j]*b[0][1];
|
||||
a[1] = b[1][2] + x1v[k]*b[1][0] + x2v[j]*b[1][1];
|
||||
a[2] = b[2][2] + x1v[k]*b[2][0] + x2v[j]*b[2][1];
|
||||
double val=a[0]*a[0]+a[1]*a[1]+a[2]*a[2];
|
||||
if (val<min)
|
||||
{
|
||||
min = val;
|
||||
x1 = x1v[k];
|
||||
x2 = x2v[j];
|
||||
}
|
||||
}
|
||||
if (x1 || x2)
|
||||
{
|
||||
b[0][2] += x1*b[0][0] + x2*b[0][1];
|
||||
b[1][2] += x1*b[1][0] + x2*b[1][1];
|
||||
b[2][2] += x1*b[2][0] + x2*b[2][1];
|
||||
r[0][2] += x1*r[0][0] + x2*r[0][1];
|
||||
r[1][2] += x1*r[1][0] + x2*r[1][1];
|
||||
r[2][2] += x1*r[2][0] + x2*r[2][1];
|
||||
greedy_recurse(b,r); // note the recursion step is here
|
||||
}
|
||||
}
|
||||
|
||||
// the meat of the greedy reduction algorithm
|
||||
void greedy_recurse(double b[3][3], int r[3][3])
|
||||
{
|
||||
col_sort(b,r);
|
||||
red12(b,r);
|
||||
red3(b,r); // recursive caller
|
||||
}
|
||||
|
||||
// set r (change of basis) to be identity then reduce basis and make it unique
|
||||
void greedy(double b[3][3],int r[3][3])
|
||||
{
|
||||
r[0][1]=r[0][2]=r[1][0]=r[1][2]=r[2][0]=r[2][1]=0;
|
||||
r[0][0]=r[1][1]=r[2][2]=1;
|
||||
greedy_recurse(b,r);
|
||||
make_unique(b,r);
|
||||
}
|
||||
|
||||
// A reduced basis isn't unique. This procedure will make it
|
||||
// "more" unique. Degenerate cases are possible, but unlikely
|
||||
// with floating point math.
|
||||
void make_unique(double b[3][3], int r[3][3])
|
||||
{
|
||||
if (fabs(b[0][0]) < fabs(b[0][1]))
|
||||
{ col_swap(b,0,1); col_swap(r,0,1); }
|
||||
if (fabs(b[0][0]) < fabs(b[0][2]))
|
||||
{ col_swap(b,0,2); col_swap(r,0,2); }
|
||||
if (fabs(b[1][1]) < fabs(b[1][2]))
|
||||
{ col_swap(b,1,2); col_swap(r,1,2); }
|
||||
|
||||
if (b[0][0] < 0){ neg_col(b,0); neg_col(r,0); }
|
||||
if (b[1][1] < 0){ neg_col(b,1); neg_col(r,1); }
|
||||
if (det(b) < 0){ neg_col(b,2); neg_col(r,2); }
|
||||
}
|
||||
}}
|
||||
131
src/USER-UEF/uef_utils.h
Normal file
131
src/USER-UEF/uef_utils.h
Normal file
@ -0,0 +1,131 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
|
||||
Contributing author: David Nicholson (MIT)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifndef LMP_UEF_UTILS_H
|
||||
#define LMP_UEF_UTILS_H
|
||||
|
||||
namespace LAMMPS_NS{ namespace UEF_utils {
|
||||
|
||||
class UEFBox
|
||||
{
|
||||
public:
|
||||
UEFBox();
|
||||
void set_strain(const double, const double);
|
||||
void step_deform(const double,const double);
|
||||
bool reduce();
|
||||
void get_box(double[3][3], double);
|
||||
void get_rot(double[3][3]);
|
||||
private:
|
||||
double l0[3][3]; // initial basis
|
||||
double w1[3],w2[3], winv[3][3]; // omega1 and omega2 (spectra of automorphisms)
|
||||
//double edot[3], delta[2];
|
||||
double theta[2];
|
||||
double l[3][3], rot[3][3], lrot[3][3];
|
||||
int r[3][3],a1[3][3],a2[3][3],a1i[3][3],a2i[3][3];
|
||||
};
|
||||
|
||||
|
||||
// lattice reduction routines
|
||||
void greedy(double[3][3],int[3][3]);
|
||||
void col_sort(double[3][3],int[3][3]);
|
||||
void red12(double[3][3],int[3][3]);
|
||||
void greedy_recurse(double[3][3],int[3][3]);
|
||||
void red3(double [3][3],int r[3][3]);
|
||||
void make_unique(double[3][3],int[3][3]);
|
||||
void rotation_matrix(double[3][3],double[3][3],const double [3][3]);
|
||||
|
||||
// A few utility functions for 3x3 arrays
|
||||
template<typename T>
|
||||
T col_prod(T x[3][3], int c1, int c2)
|
||||
{
|
||||
return x[0][c1]*x[0][c2]+x[1][c1]*x[1][c2]+x[2][c1]*x[2][c2];
|
||||
}
|
||||
|
||||
template<typename T>
|
||||
void col_swap(T x[3][3], int c1, int c2)
|
||||
{
|
||||
for (int k=0;k<3;k++)
|
||||
{
|
||||
T t = x[k][c2];
|
||||
x[k][c2]=x[k][c1];
|
||||
x[k][c1]=t;
|
||||
}
|
||||
}
|
||||
|
||||
template<typename T>
|
||||
void neg_col(T x[3][3], int c1)
|
||||
{
|
||||
x[0][c1] = -x[0][c1];
|
||||
x[1][c1] = -x[1][c1];
|
||||
x[2][c1] = -x[2][c1];
|
||||
}
|
||||
|
||||
template<typename T>
|
||||
void neg_row(T x[3][3], int c1)
|
||||
{
|
||||
x[c1][0] = -x[c1][0];
|
||||
x[c1][1] = -x[c1][1];
|
||||
x[c1][2] = -x[c1][2];
|
||||
}
|
||||
|
||||
template<typename T>
|
||||
T det(T x[3][3])
|
||||
{
|
||||
double val;
|
||||
val = x[0][0]*(x[1][1]*x[2][2] - x[1][2]*x[2][1]);
|
||||
val -= x[0][1]*(x[1][0]*x[2][2] - x[1][2]*x[2][0]);
|
||||
val += x[0][2]*(x[1][0]*x[2][1] - x[1][1]*x[2][0]);
|
||||
return val;
|
||||
}
|
||||
|
||||
template<typename T>
|
||||
bool mat_same(T x1[3][3], T x2[3][3])
|
||||
{
|
||||
bool v = true;
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
v &= (x1[k][j]==x2[k][j]);
|
||||
return v;
|
||||
}
|
||||
|
||||
template<typename T>
|
||||
void mul_m1(T m1[3][3], const T m2[3][3])
|
||||
{
|
||||
T t[3][3];
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
t[k][j]=m1[k][j];
|
||||
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
m1[k][j] = t[k][0]*m2[0][j] + t[k][1]*m2[1][j] + t[k][2]*m2[2][j];
|
||||
}
|
||||
|
||||
template<typename T>
|
||||
void mul_m2(const T m1[3][3], T m2[3][3])
|
||||
{
|
||||
T t[3][3];
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
t[k][j]=m2[k][j];
|
||||
|
||||
for (int k=0;k<3;k++)
|
||||
for (int j=0;j<3;j++)
|
||||
m2[k][j] = m1[k][0]*t[0][j] + m1[k][1]*t[1][j] + m1[k][2]*t[2][j];
|
||||
}
|
||||
|
||||
}
|
||||
}
|
||||
#endif
|
||||
@ -353,6 +353,15 @@ FixNH::FixNH(LAMMPS *lmp, int narg, char **arg) :
|
||||
} else if (strcmp(arg[iarg],"disc") == 0) {
|
||||
iarg++;
|
||||
|
||||
// keywords erate, strain, and ext are also parsed in fix/nh/uef
|
||||
|
||||
} else if (strcmp(arg[iarg],"erate") == 0) {
|
||||
iarg += 3;
|
||||
} else if (strcmp(arg[iarg],"strain") == 0) {
|
||||
iarg += 3;
|
||||
} else if (strcmp(arg[iarg],"ext") == 0) {
|
||||
iarg += 2;
|
||||
|
||||
} else error->all(FLERR,"Illegal fix nvt/npt/nph command");
|
||||
}
|
||||
|
||||
|
||||
Reference in New Issue
Block a user