LAMMPS Users Manual

diff --git a/doc/src/Manual.txt b/doc/src/Manual.txt index 25319570a..3e2c7dc82 100644 --- a/doc/src/Manual.txt +++ b/doc/src/Manual.txt @@ -1,340 +1,340 @@ LAMMPS Users Manual - + "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

LAMMPS Documentation :c,h3 -26 Jan 2017 version :c,h4 +10 Feb 2016 version :c,h4 Version info: :h4 The LAMMPS "version" is the date when it was released, such as 1 May 2010. LAMMPS is updated continuously. Whenever we fix a bug or add a feature, we release it immediately, and post a notice on "this page of the WWW site"_bug. Every 2-4 months one of the incremental releases is subjected to more thorough testing and labeled as a {stable} version. Each dated copy of LAMMPS contains all the features and bug-fixes up to and including that version date. The version date is printed to the screen and logfile every time you run LAMMPS. It is also in the file src/version.h and in the LAMMPS directory name created when you unpack a tarball, and at the top of the first page of the manual (this page). If you browse the HTML doc pages on the LAMMPS WWW site, they always describe the most current version of LAMMPS. :ulb,l If you browse the HTML doc pages included in your tarball, they describe the version you have. :l The "PDF file"_Manual.pdf on the WWW site or in the tarball is updated about once per month. This is because it is large, and we don't want it to be part of every patch. :l There is also a "Developer.pdf"_Developer.pdf file in the doc directory, which describes the internal structure and algorithms of LAMMPS. :l :ule LAMMPS stands for Large-scale Atomic/Molecular Massively Parallel Simulator. LAMMPS is a classical molecular dynamics simulation code designed to run efficiently on parallel computers. It was developed at Sandia National Laboratories, a US Department of Energy facility, with funding from the DOE. It is an open-source code, distributed freely under the terms of the GNU Public License (GPL). The current core group of LAMMPS developers is at Sandia National Labs and Temple University: "Steve Plimpton"_sjp, sjplimp at sandia.gov :ulb,l Aidan Thompson, athomps at sandia.gov :l Stan Moore, stamoore at sandia.gov :l "Axel Kohlmeyer"_ako, akohlmey at gmail.com :l :ule Past core developers include Paul Crozier, Ray Shan and Mark Stevens, all at Sandia. The [LAMMPS home page] at "http://lammps.sandia.gov"_http://lammps.sandia.gov has more information about the code and its uses. Interaction with external LAMMPS developers, bug reports and feature requests are mainly coordinated through the "LAMMPS project on GitHub."_https://github.com/lammps/lammps The lammps.org domain, currently hosting "public continuous integration testing"_https://ci.lammps.org/job/lammps/ and "precompiled Linux RPM and Windows installer packages"_http://rpm.lammps.org is located at Temple University and managed by Richard Berger, richard.berger at temple.edu. :link(bug,http://lammps.sandia.gov/bug.html) :link(sjp,http://www.sandia.gov/~sjplimp) :link(ako,http://goo.gl/1wk0) :line The LAMMPS documentation is organized into the following sections. If you find errors or omissions in this manual or have suggestions for useful information to add, please send an email to the developers so we can improve the LAMMPS documentation. Once you are familiar with LAMMPS, you may want to bookmark "this page"_Section_commands.html#comm at Section_commands.html#comm since it gives quick access to documentation for all LAMMPS commands. "PDF file"_Manual.pdf of the entire manual, generated by "htmldoc"_http://freecode.com/projects/htmldoc "Introduction"_Section_intro.html :olb,l 1.1 "What is LAMMPS"_intro_1 :ulb,b 1.2 "LAMMPS features"_intro_2 :b 1.3 "LAMMPS non-features"_intro_3 :b 1.4 "Open source distribution"_intro_4 :b 1.5 "Acknowledgments and citations"_intro_5 :ule,b "Getting started"_Section_start.html :l 2.1 "What's in the LAMMPS distribution"_start_1 :ulb,b 2.2 "Making LAMMPS"_start_2 :b 2.3 "Making LAMMPS with optional packages"_start_3 :b 2.4 "Building LAMMPS via the Make.py script"_start_4 :b 2.5 "Building LAMMPS as a library"_start_5 :b 2.6 "Running LAMMPS"_start_6 :b 2.7 "Command-line options"_start_7 :b 2.8 "Screen output"_start_8 :b 2.9 "Tips for users of previous versions"_start_9 :ule,b "Commands"_Section_commands.html :l 3.1 "LAMMPS input script"_cmd_1 :ulb,b 3.2 "Parsing rules"_cmd_2 :b 3.3 "Input script structure"_cmd_3 :b 3.4 "Commands listed by category"_cmd_4 :b 3.5 "Commands listed alphabetically"_cmd_5 :ule,b "Packages"_Section_packages.html :l 4.1 "Standard packages"_pkg_1 :ulb,b 4.2 "User packages"_pkg_2 :ule,b "Accelerating LAMMPS performance"_Section_accelerate.html :l 5.1 "Measuring performance"_acc_1 :ulb,b 5.2 "Algorithms and code options to boost performace"_acc_2 :b 5.3 "Accelerator packages with optimized styles"_acc_3 :b 5.3.1 "GPU package"_accelerate_gpu.html :ulb,b 5.3.2 "USER-INTEL package"_accelerate_intel.html :b 5.3.3 "KOKKOS package"_accelerate_kokkos.html :b 5.3.4 "USER-OMP package"_accelerate_omp.html :b 5.3.5 "OPT package"_accelerate_opt.html :ule,b 5.4 "Comparison of various accelerator packages"_acc_4 :ule,b "How-to discussions"_Section_howto.html :l 6.1 "Restarting a simulation"_howto_1 :ulb,b 6.2 "2d simulations"_howto_2 :b 6.3 "CHARMM and AMBER force fields"_howto_3 :b 6.4 "Running multiple simulations from one input script"_howto_4 :b 6.5 "Multi-replica simulations"_howto_5 :b 6.6 "Granular models"_howto_6 :b 6.7 "TIP3P water model"_howto_7 :b 6.8 "TIP4P water model"_howto_8 :b 6.9 "SPC water model"_howto_9 :b 6.10 "Coupling LAMMPS to other codes"_howto_10 :b 6.11 "Visualizing LAMMPS snapshots"_howto_11 :b 6.12 "Triclinic (non-orthogonal) simulation boxes"_howto_12 :b 6.13 "NEMD simulations"_howto_13 :b 6.14 "Finite-size spherical and aspherical particles"_howto_14 :b 6.15 "Output from LAMMPS (thermo, dumps, computes, fixes, variables)"_howto_15 :b 6.16 "Thermostatting, barostatting, and compute temperature"_howto_16 :b 6.17 "Walls"_howto_17 :b 6.18 "Elastic constants"_howto_18 :b 6.19 "Library interface to LAMMPS"_howto_19 :b 6.20 "Calculating thermal conductivity"_howto_20 :b 6.21 "Calculating viscosity"_howto_21 :b 6.22 "Calculating a diffusion coefficient"_howto_22 :b 6.23 "Using chunks to calculate system properties"_howto_23 :b 6.24 "Setting parameters for pppm/disp"_howto_24 :b 6.25 "Polarizable models"_howto_25 :b 6.26 "Adiabatic core/shell model"_howto_26 :b 6.27 "Drude induced dipoles"_howto_27 :ule,b "Example problems"_Section_example.html :l "Performance & scalability"_Section_perf.html :l "Additional tools"_Section_tools.html :l "Modifying & extending LAMMPS"_Section_modify.html :l 10.1 "Atom styles"_mod_1 :ulb,b 10.2 "Bond, angle, dihedral, improper potentials"_mod_2 :b 10.3 "Compute styles"_mod_3 :b 10.4 "Dump styles"_mod_4 :b 10.5 "Dump custom output options"_mod_5 :b 10.6 "Fix styles"_mod_6 :b 10.7 "Input script commands"_mod_7 :b 10.8 "Kspace computations"_mod_8 :b 10.9 "Minimization styles"_mod_9 :b 10.10 "Pairwise potentials"_mod_10 :b 10.11 "Region styles"_mod_11 :b 10.12 "Body styles"_mod_12 :b 10.13 "Thermodynamic output options"_mod_13 :b 10.14 "Variable options"_mod_14 :b 10.15 "Submitting new features for inclusion in LAMMPS"_mod_15 :ule,b "Python interface"_Section_python.html :l 11.1 "Overview of running LAMMPS from Python"_py_1 :ulb,b 11.2 "Overview of using Python from a LAMMPS script"_py_2 :b 11.3 "Building LAMMPS as a shared library"_py_3 :b 11.4 "Installing the Python wrapper into Python"_py_4 :b 11.5 "Extending Python with MPI to run in parallel"_py_5 :b 11.6 "Testing the Python-LAMMPS interface"_py_6 :b 11.7 "Using LAMMPS from Python"_py_7 :b 11.8 "Example Python scripts that use LAMMPS"_py_8 :ule,b "Errors"_Section_errors.html :l 12.1 "Common problems"_err_1 :ulb,b 12.2 "Reporting bugs"_err_2 :b 12.3 "Error & warning messages"_err_3 :ule,b "Future and history"_Section_history.html :l 13.1 "Coming attractions"_hist_1 :ulb,b 13.2 "Past versions"_hist_2 :ule,b :ole :link(intro_1,Section_intro.html#intro_1) :link(intro_2,Section_intro.html#intro_2) :link(intro_3,Section_intro.html#intro_3) :link(intro_4,Section_intro.html#intro_4) :link(intro_5,Section_intro.html#intro_5) :link(start_1,Section_start.html#start_1) :link(start_2,Section_start.html#start_2) :link(start_3,Section_start.html#start_3) :link(start_4,Section_start.html#start_4) :link(start_5,Section_start.html#start_5) :link(start_6,Section_start.html#start_6) :link(start_7,Section_start.html#start_7) :link(start_8,Section_start.html#start_8) :link(start_9,Section_start.html#start_9) :link(cmd_1,Section_commands.html#cmd_1) :link(cmd_2,Section_commands.html#cmd_2) :link(cmd_3,Section_commands.html#cmd_3) :link(cmd_4,Section_commands.html#cmd_4) :link(cmd_5,Section_commands.html#cmd_5) :link(pkg_1,Section_packages.html#pkg_1) :link(pkg_2,Section_packages.html#pkg_2) :link(acc_1,Section_accelerate.html#acc_1) :link(acc_2,Section_accelerate.html#acc_2) :link(acc_3,Section_accelerate.html#acc_3) :link(acc_4,Section_accelerate.html#acc_4) :link(howto_1,Section_howto.html#howto_1) :link(howto_2,Section_howto.html#howto_2) :link(howto_3,Section_howto.html#howto_3) :link(howto_4,Section_howto.html#howto_4) :link(howto_5,Section_howto.html#howto_5) :link(howto_6,Section_howto.html#howto_6) :link(howto_7,Section_howto.html#howto_7) :link(howto_8,Section_howto.html#howto_8) :link(howto_9,Section_howto.html#howto_9) :link(howto_10,Section_howto.html#howto_10) :link(howto_11,Section_howto.html#howto_11) :link(howto_12,Section_howto.html#howto_12) :link(howto_13,Section_howto.html#howto_13) :link(howto_14,Section_howto.html#howto_14) :link(howto_15,Section_howto.html#howto_15) :link(howto_16,Section_howto.html#howto_16) :link(howto_17,Section_howto.html#howto_17) :link(howto_18,Section_howto.html#howto_18) :link(howto_19,Section_howto.html#howto_19) :link(howto_20,Section_howto.html#howto_20) :link(howto_21,Section_howto.html#howto_21) :link(howto_22,Section_howto.html#howto_22) :link(howto_23,Section_howto.html#howto_23) :link(howto_24,Section_howto.html#howto_24) :link(howto_25,Section_howto.html#howto_25) :link(howto_26,Section_howto.html#howto_26) :link(howto_27,Section_howto.html#howto_27) :link(mod_1,Section_modify.html#mod_1) :link(mod_2,Section_modify.html#mod_2) :link(mod_3,Section_modify.html#mod_3) :link(mod_4,Section_modify.html#mod_4) :link(mod_5,Section_modify.html#mod_5) :link(mod_6,Section_modify.html#mod_6) :link(mod_7,Section_modify.html#mod_7) :link(mod_8,Section_modify.html#mod_8) :link(mod_9,Section_modify.html#mod_9) :link(mod_10,Section_modify.html#mod_10) :link(mod_11,Section_modify.html#mod_11) :link(mod_12,Section_modify.html#mod_12) :link(mod_13,Section_modify.html#mod_13) :link(mod_14,Section_modify.html#mod_14) :link(mod_15,Section_modify.html#mod_15) :link(py_1,Section_python.html#py_1) :link(py_2,Section_python.html#py_2) :link(py_3,Section_python.html#py_3) :link(py_4,Section_python.html#py_4) :link(py_5,Section_python.html#py_5) :link(py_6,Section_python.html#py_6) :link(err_1,Section_errors.html#err_1) :link(err_2,Section_errors.html#err_2) :link(err_3,Section_errors.html#err_3) :link(hist_1,Section_history.html#hist_1) :link(hist_2,Section_history.html#hist_2) diff --git a/doc/src/Section_commands.txt b/doc/src/Section_commands.txt index 7284b4b0f..30c2743ab 100644 --- a/doc/src/Section_commands.txt +++ b/doc/src/Section_commands.txt @@ -1,1215 +1,1215 @@ "Previous Section"_Section_start.html - "LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc - "Next Section"_Section_packages.html :c :link(lws,http://lammps.sandia.gov) :link(ld,Manual.html) :link(lc,Section_commands.html#comm) :line 3. Commands :h3 This section describes how a LAMMPS input script is formatted and the input script commands used to define a LAMMPS simulation. 3.1 "LAMMPS input script"_#cmd_1 3.2 "Parsing rules"_#cmd_2 3.3 "Input script structure"_#cmd_3 3.4 "Commands listed by category"_#cmd_4 3.5 "Commands listed alphabetically"_#cmd_5 :all(b) :line :line 3.1 LAMMPS input script :link(cmd_1),h4 LAMMPS executes by reading commands from a input script (text file), one line at a time. When the input script ends, LAMMPS exits. Each command causes LAMMPS to take some action. It may set an internal variable, read in a file, or run a simulation. Most commands have default settings, which means you only need to use the command if you wish to change the default. In many cases, the ordering of commands in an input script is not important. However the following rules apply: (1) LAMMPS does not read your entire input script and then perform a simulation with all the settings. Rather, the input script is read one line at a time and each command takes effect when it is read. Thus this sequence of commands: timestep 0.5 run 100 run 100 :pre does something different than this sequence: run 100 timestep 0.5 run 100 :pre In the first case, the specified timestep (0.5 fmsec) is used for two simulations of 100 timesteps each. In the 2nd case, the default timestep (1.0 fmsec) is used for the 1st 100 step simulation and a 0.5 fmsec timestep is used for the 2nd one. (2) Some commands are only valid when they follow other commands. For example you cannot set the temperature of a group of atoms until atoms have been defined and a group command is used to define which atoms belong to the group. (3) Sometimes command B will use values that can be set by command A. This means command A must precede command B in the input script if it is to have the desired effect. For example, the "read_data"_read_data.html command initializes the system by setting up the simulation box and assigning atoms to processors. If default values are not desired, the "processors"_processors.html and "boundary"_boundary.html commands need to be used before read_data to tell LAMMPS how to map processors to the simulation box. Many input script errors are detected by LAMMPS and an ERROR or WARNING message is printed. "This section"_Section_errors.html gives more information on what errors mean. The documentation for each command lists restrictions on how the command can be used. :line 3.2 Parsing rules :link(cmd_2),h4 Each non-blank line in the input script is treated as a command. LAMMPS commands are case sensitive. Command names are lower-case, as are specified command arguments. Upper case letters may be used in file names or user-chosen ID strings. Here is how each line in the input script is parsed by LAMMPS: (1) If the last printable character on the line is a "&" character, the command is assumed to continue on the next line. The next line is concatenated to the previous line by removing the "&" character and line break. This allows long commands to be continued across two or more lines. See the discussion of triple quotes in (6) for how to continue a command across multiple line without using "&" characters. (2) All characters from the first "#" character onward are treated as comment and discarded. See an exception in (6). Note that a comment after a trailing "&" character will prevent the command from continuing on the next line. Also note that for multi-line commands a single leading "#" will comment out the entire command. (3) The line is searched repeatedly for $ characters, which indicate variables that are replaced with a text string. See an exception in (6). If the $ is followed by curly brackets, then the variable name is the text inside the curly brackets. If no curly brackets follow the $, then the variable name is the single character immediately following the $. Thus $\{myTemp\} and $x refer to variable names "myTemp" and "x". How the variable is converted to a text string depends on what style of variable it is; see the "variable"_variable.html doc page for details. It can be a variable that stores multiple text strings, and return one of them. The returned text string can be multiple "words" (space separated) which will then be interpreted as multiple arguments in the input command. The variable can also store a numeric formula which will be evaluated and its numeric result returned as a string. As a special case, if the $ is followed by parenthesis, then the text inside the parenthesis is treated as an "immediate" variable and evaluated as an "equal-style variable"_variable.html. This is a way to use numeric formulas in an input script without having to assign them to variable names. For example, these 3 input script lines: variable X equal (xlo+xhi)/2+sqrt(v_area) region 1 block $X 2 INF INF EDGE EDGE variable X delete :pre can be replaced by region 1 block $((xlo+xhi)/2+sqrt(v_area)) 2 INF INF EDGE EDGE :pre so that you do not have to define (or discard) a temporary variable X. Note that neither the curly-bracket or immediate form of variables can contain nested $ characters for other variables to substitute for. Thus you cannot do this: variable a equal 2 variable b2 equal 4 print "B2 = $\{b$a\}" :pre Nor can you specify this $($x-1.0) for an immediate variable, but you could use $(v_x-1.0), since the latter is valid syntax for an "equal-style variable"_variable.html. See the "variable"_variable.html command for more details of how strings are assigned to variables and evaluated, and how they can be used in input script commands. (4) The line is broken into "words" separated by whitespace (tabs, spaces). Note that words can thus contain letters, digits, underscores, or punctuation characters. (5) The first word is the command name. All successive words in the line are arguments. (6) If you want text with spaces to be treated as a single argument, it can be enclosed in either single or double or triple quotes. A long single argument enclosed in single or double quotes can span multiple lines if the "&" character is used, as described above. When the lines are concatenated together (and the "&" characters and line breaks removed), the text will become a single line. If you want multiple lines of an argument to retain their line breaks, the text can be enclosed in triple quotes, in which case "&" characters are not needed. For example: print "Volume = $v" print 'Volume = $v' if "$\{steps\} > 1000" then quit variable a string "red green blue & purple orange cyan" print """ System volume = $v System temperature = $t """ :pre In each case, the single, double, or triple quotes are removed when the single argument they enclose is stored internally. See the "dump modify format"_dump_modify.html, "print"_print.html, "if"_if.html, and "python"_python.html commands for examples. A "#" or "$" character that is between quotes will not be treated as a comment indicator in (2) or substituted for as a variable in (3). NOTE: If the argument is itself a command that requires a quoted argument (e.g. using a "print"_print.html command as part of an "if"_if.html or "run every"_run.html command), then single, double, or triple quotes can be nested in the usual manner. See the doc pages for those commands for examples. Only one of level of nesting is allowed, but that should be sufficient for most use cases. :line 3.3 Input script structure :h4,link(cmd_3) This section describes the structure of a typical LAMMPS input script. The "examples" directory in the LAMMPS distribution contains many sample input scripts; the corresponding problems are discussed in "Section 7"_Section_example.html, and animated on the "LAMMPS WWW Site"_lws. A LAMMPS input script typically has 4 parts: Initialization Atom definition Settings Run a simulation :ol The last 2 parts can be repeated as many times as desired. I.e. run a simulation, change some settings, run some more, etc. Each of the 4 parts is now described in more detail. Remember that almost all the commands need only be used if a non-default value is desired. (1) Initialization Set parameters that need to be defined before atoms are created or read-in from a file. The relevant commands are "units"_units.html, "dimension"_dimension.html, "newton"_newton.html, "processors"_processors.html, "boundary"_boundary.html, "atom_style"_atom_style.html, "atom_modify"_atom_modify.html. If force-field parameters appear in the files that will be read, these commands tell LAMMPS what kinds of force fields are being used: "pair_style"_pair_style.html, "bond_style"_bond_style.html, "angle_style"_angle_style.html, "dihedral_style"_dihedral_style.html, "improper_style"_improper_style.html. (2) Atom definition There are 3 ways to define atoms in LAMMPS. Read them in from a data or restart file via the "read_data"_read_data.html or "read_restart"_read_restart.html commands. These files can contain molecular topology information. Or create atoms on a lattice (with no molecular topology), using these commands: "lattice"_lattice.html, "region"_region.html, "create_box"_create_box.html, "create_atoms"_create_atoms.html. The entire set of atoms can be duplicated to make a larger simulation using the "replicate"_replicate.html command. (3) Settings Once atoms and molecular topology are defined, a variety of settings can be specified: force field coefficients, simulation parameters, output options, etc. Force field coefficients are set by these commands (they can also be set in the read-in files): "pair_coeff"_pair_coeff.html, "bond_coeff"_bond_coeff.html, "angle_coeff"_angle_coeff.html, "dihedral_coeff"_dihedral_coeff.html, "improper_coeff"_improper_coeff.html, "kspace_style"_kspace_style.html, "dielectric"_dielectric.html, "special_bonds"_special_bonds.html. Various simulation parameters are set by these commands: "neighbor"_neighbor.html, "neigh_modify"_neigh_modify.html, "group"_group.html, "timestep"_timestep.html, "reset_timestep"_reset_timestep.html, "run_style"_run_style.html, "min_style"_min_style.html, "min_modify"_min_modify.html. Fixes impose a variety of boundary conditions, time integration, and diagnostic options. The "fix"_fix.html command comes in many flavors. Various computations can be specified for execution during a simulation using the "compute"_compute.html, "compute_modify"_compute_modify.html, and "variable"_variable.html commands. Output options are set by the "thermo"_thermo.html, "dump"_dump.html, and "restart"_restart.html commands. (4) Run a simulation A molecular dynamics simulation is run using the "run"_run.html command. Energy minimization (molecular statics) is performed using the "minimize"_minimize.html command. A parallel tempering (replica-exchange) simulation can be run using the "temper"_temper.html command. :line 3.4 Commands listed by category :link(cmd_4),h4 This section lists all LAMMPS commands, grouped by category. The "next section"_#cmd_5 lists the same commands alphabetically. The next section also includes (long) lists of style options for entries that appear in the following categories as a single command (fix, compute, pair, etc). Commands that are added by user packages are not included in these categories, but they are in the next section. Initialization: "newton"_newton.html, "package"_package.html, "processors"_processors.html, "suffix"_suffix.html, "units"_units.html Setup simulation box: "boundary"_boundary.html, "box"_box.html, "change_box"_change_box.html, "create_box"_create_box.html, "dimension"_dimension.html, "lattice"_lattice.html, "region"_region.html Setup atoms: "atom_modify"_atom_modify.html, "atom_style"_atom_style.html, "balance"_balance.html, "create_atoms"_create_atoms.html, "create_bonds"_create_bonds.html, "delete_atoms"_delete_atoms.html, "delete_bonds"_delete_bonds.html, "displace_atoms"_displace_atoms.html, "group"_group.html, "mass"_mass.html, "molecule"_molecule.html, "read_data"_read_data.html, "read_dump"_read_dump.html, "read_restart"_read_restart.html, "replicate"_replicate.html, "set"_set.html, "velocity"_velocity.html Force fields: "angle_coeff"_angle_coeff.html, "angle_style"_angle_style.html, "bond_coeff"_bond_coeff.html, "bond_style"_bond_style.html, "bond_write"_bond_write.html, "dielectric"_dielectric.html, "dihedral_coeff"_dihedral_coeff.html, "dihedral_style"_dihedral_style.html, "improper_coeff"_improper_coeff.html, "improper_style"_improper_style.html, "kspace_modify"_kspace_modify.html, "kspace_style"_kspace_style.html, "pair_coeff"_pair_coeff.html, "pair_modify"_pair_modify.html, "pair_style"_pair_style.html, "pair_write"_pair_write.html, "special_bonds"_special_bonds.html Settings: "comm_modify"_comm_modify.html, "comm_style"_comm_style.html, "info"_info.html, "min_modify"_min_modify.html, "min_style"_min_style.html, "neigh_modify"_neigh_modify.html, "neighbor"_neighbor.html, "partition"_partition.html, "reset_timestep"_reset_timestep.html, "run_style"_run_style.html, "timer"_timer.html, "timestep"_timestep.html Operations within timestepping (fixes) and diagnositics (computes): "compute"_compute.html, "compute_modify"_compute_modify.html, "fix"_fix.html, "fix_modify"_fix_modify.html, "uncompute"_uncompute.html, "unfix"_unfix.html Output: "dump image"_dump_image.html, "dump movie"_dump_image.html, "dump"_dump.html, "dump_modify"_dump_modify.html, "restart"_restart.html, "thermo"_thermo.html, "thermo_modify"_thermo_modify.html, "thermo_style"_thermo_style.html, "undump"_undump.html, "write_coeff"_write_coeff.html, "write_data"_write_data.html, "write_dump"_write_dump.html, "write_restart"_write_restart.html Actions: "minimize"_minimize.html, "neb"_neb.html, "prd"_prd.html, "rerun"_rerun.html, "run"_run.html, "tad"_tad.html, "temper"_temper.html Input script control: "clear"_clear.html, "echo"_echo.html, "if"_if.html, "include"_include.html, "jump"_jump.html, "label"_label.html, "log"_log.html, "next"_next.html, "print"_print.html, "python"_python.html, "quit"_quit.html, "shell"_shell.html, "variable"_variable.html :line 3.5 Individual commands :h4,link(cmd_5),link(comm) This section lists all LAMMPS commands alphabetically, with a separate listing below of styles within certain commands. The "previous section"_#cmd_4 lists the same commands, grouped by category. Note that some style options for some commands are part of specific LAMMPS packages, which means they cannot be used unless the package was included when LAMMPS was built. Not all packages are included in a default LAMMPS build. These dependencies are listed as Restrictions in the command's documentation. "angle_coeff"_angle_coeff.html, "angle_style"_angle_style.html, "atom_modify"_atom_modify.html, "atom_style"_atom_style.html, "balance"_balance.html, "bond_coeff"_bond_coeff.html, "bond_style"_bond_style.html, "bond_write"_bond_write.html, "boundary"_boundary.html, "box"_box.html, "change_box"_change_box.html, "clear"_clear.html, "comm_modify"_comm_modify.html, "comm_style"_comm_style.html, "compute"_compute.html, "compute_modify"_compute_modify.html, "create_atoms"_create_atoms.html, "create_bonds"_create_bonds.html, "create_box"_create_box.html, "delete_atoms"_delete_atoms.html, "delete_bonds"_delete_bonds.html, "dielectric"_dielectric.html, "dihedral_coeff"_dihedral_coeff.html, "dihedral_style"_dihedral_style.html, "dimension"_dimension.html, "displace_atoms"_displace_atoms.html, "dump"_dump.html, "dump image"_dump_image.html, "dump_modify"_dump_modify.html, "dump movie"_dump_image.html, "echo"_echo.html, "fix"_fix.html, "fix_modify"_fix_modify.html, "group"_group.html, "if"_if.html, "info"_info.html, "improper_coeff"_improper_coeff.html, "improper_style"_improper_style.html, "include"_include.html, "jump"_jump.html, "kspace_modify"_kspace_modify.html, "kspace_style"_kspace_style.html, "label"_label.html, "lattice"_lattice.html, "log"_log.html, "mass"_mass.html, "minimize"_minimize.html, "min_modify"_min_modify.html, "min_style"_min_style.html, "molecule"_molecule.html, "neb"_neb.html, "neigh_modify"_neigh_modify.html, "neighbor"_neighbor.html, "newton"_newton.html, "next"_next.html, "package"_package.html, "pair_coeff"_pair_coeff.html, "pair_modify"_pair_modify.html, "pair_style"_pair_style.html, "pair_write"_pair_write.html, "partition"_partition.html, "prd"_prd.html, "print"_print.html, "processors"_processors.html, "python"_python.html, "quit"_quit.html, "read_data"_read_data.html, "read_dump"_read_dump.html, "read_restart"_read_restart.html, "region"_region.html, "replicate"_replicate.html, "rerun"_rerun.html, "reset_timestep"_reset_timestep.html, "restart"_restart.html, "run"_run.html, "run_style"_run_style.html, "set"_set.html, "shell"_shell.html, "special_bonds"_special_bonds.html, "suffix"_suffix.html, "tad"_tad.html, "temper"_temper.html, "thermo"_thermo.html, "thermo_modify"_thermo_modify.html, "thermo_style"_thermo_style.html, "timer"_timer.html, "timestep"_timestep.html, "uncompute"_uncompute.html, "undump"_undump.html, "unfix"_unfix.html, "units"_units.html, "variable"_variable.html, "velocity"_velocity.html, "write_coeff"_write_coeff.html, "write_data"_write_data.html, "write_dump"_write_dump.html, "write_restart"_write_restart.html :tb(c=6,ea=c) These are additional commands in USER packages, which can be used if "LAMMPS is built with the appropriate package"_Section_start.html#start_3. "dump custom/vtk"_dump_custom_vtk.html, "dump nc"_dump_nc.html, "dump nc/mpiio"_dump_nc.html, "group2ndx"_group2ndx.html, "ndx2group"_group2ndx.html, "temper/grem"_temper_grem.html :tb(c=3,ea=c) :line Fix styles :h4 See the "fix"_fix.html command for one-line descriptions of each style or click on the style itself for a full description. Some of the styles have accelerated versions, which can be used if LAMMPS is built with the "appropriate accelerated package"_Section_accelerate.html. This is indicated by additional letters in parenthesis: g = GPU, i = USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT. "adapt"_fix_adapt.html, "addforce"_fix_addforce.html, "append/atoms"_fix_append_atoms.html, "atom/swap"_fix_atom_swap.html, "aveforce"_fix_aveforce.html, "ave/atom"_fix_ave_atom.html, "ave/chunk"_fix_ave_chunk.html, "ave/correlate"_fix_ave_correlate.html, "ave/histo"_fix_ave_histo.html, "ave/histo/weight"_fix_ave_histo.html, "ave/time"_fix_ave_time.html, "balance"_fix_balance.html, "bond/break"_fix_bond_break.html, "bond/create"_fix_bond_create.html, "bond/swap"_fix_bond_swap.html, "box/relax"_fix_box_relax.html, "cmap"_fix_cmap.html, "controller"_fix_controller.html, "deform (k)"_fix_deform.html, "deposit"_fix_deposit.html, "drag"_fix_drag.html, "dt/reset"_fix_dt_reset.html, "efield"_fix_efield.html, "ehex"_fix_ehex.html, "enforce2d"_fix_enforce2d.html, "evaporate"_fix_evaporate.html, "external"_fix_external.html, "freeze"_fix_freeze.html, "gcmc"_fix_gcmc.html, "gld"_fix_gld.html, "gravity (o)"_fix_gravity.html, "halt"_fix_halt.html, "heat"_fix_heat.html, "indent"_fix_indent.html, "langevin (k)"_fix_langevin.html, "lineforce"_fix_lineforce.html, "momentum (k)"_fix_momentum.html, "move"_fix_move.html, "mscg"_fix_mscg.html, "msst"_fix_msst.html, "neb"_fix_neb.html, "nph (ko)"_fix_nh.html, "nphug (o)"_fix_nphug.html, "nph/asphere (o)"_fix_nph_asphere.html, "nph/body"_fix_nph_body.html, "nph/sphere (o)"_fix_nph_sphere.html, "npt (kio)"_fix_nh.html, "npt/asphere (o)"_fix_npt_asphere.html, "npt/body"_fix_npt_body.html, "npt/sphere (o)"_fix_npt_sphere.html, "nve (kio)"_fix_nve.html, "nve/asphere (i)"_fix_nve_asphere.html, "nve/asphere/noforce"_fix_nve_asphere_noforce.html, "nve/body"_fix_nve_body.html, "nve/limit"_fix_nve_limit.html, "nve/line"_fix_nve_line.html, "nve/noforce"_fix_nve_noforce.html, "nve/sphere (o)"_fix_nve_sphere.html, "nve/tri"_fix_nve_tri.html, "nvt (iko)"_fix_nh.html, "nvt/asphere (o)"_fix_nvt_asphere.html, "nvt/body"_fix_nvt_body.html, "nvt/sllod (io)"_fix_nvt_sllod.html, "nvt/sphere (o)"_fix_nvt_sphere.html, "oneway"_fix_oneway.html, "orient/bcc"_fix_orient.html, "orient/fcc"_fix_orient.html, "planeforce"_fix_planeforce.html, "poems"_fix_poems.html, "pour"_fix_pour.html, "press/berendsen"_fix_press_berendsen.html, "print"_fix_print.html, "property/atom"_fix_property_atom.html, "qeq/comb (o)"_fix_qeq_comb.html, "qeq/dynamic"_fix_qeq.html, "qeq/fire"_fix_qeq.html, "qeq/point"_fix_qeq.html, "qeq/shielded"_fix_qeq.html, "qeq/slater"_fix_qeq.html, "rattle"_fix_shake.html, "reax/bonds"_fix_reax_bonds.html, "recenter"_fix_recenter.html, "restrain"_fix_restrain.html, "rigid (o)"_fix_rigid.html, "rigid/nph (o)"_fix_rigid.html, "rigid/npt (o)"_fix_rigid.html, "rigid/nve (o)"_fix_rigid.html, "rigid/nvt (o)"_fix_rigid.html, "rigid/small (o)"_fix_rigid.html, "rigid/small/nph (o)"_fix_rigid.html, "rigid/small/npt (o)"_fix_rigid.html, "rigid/small/nve (o)"_fix_rigid.html, "rigid/small/nvt (o)"_fix_rigid.html, "setforce (k)"_fix_setforce.html, "shake"_fix_shake.html, "spring"_fix_spring.html, "spring/chunk"_fix_spring_chunk.html, "spring/rg"_fix_spring_rg.html, "spring/self"_fix_spring_self.html, "srd"_fix_srd.html, "store/force"_fix_store_force.html, "store/state"_fix_store_state.html, "temp/berendsen"_fix_temp_berendsen.html, "temp/csld"_fix_temp_csvr.html, "temp/csvr"_fix_temp_csvr.html, "temp/rescale"_fix_temp_rescale.html, "tfmc"_fix_tfmc.html, "thermal/conductivity"_fix_thermal_conductivity.html, "tmd"_fix_tmd.html, "ttm"_fix_ttm.html, "tune/kspace"_fix_tune_kspace.html, "vector"_fix_vector.html, "viscosity"_fix_viscosity.html, "viscous"_fix_viscous.html, "wall/colloid"_fix_wall.html, "wall/gran"_fix_wall_gran.html, "wall/gran/region"_fix_wall_gran_region.html, "wall/harmonic"_fix_wall.html, "wall/lj1043"_fix_wall.html, "wall/lj126"_fix_wall.html, "wall/lj93"_fix_wall.html, "wall/piston"_fix_wall_piston.html, "wall/reflect (k)"_fix_wall_reflect.html, "wall/region"_fix_wall_region.html, "wall/srd"_fix_wall_srd.html :tb(c=8,ea=c) These are additional fix styles in USER packages, which can be used if "LAMMPS is built with the appropriate package"_Section_start.html#start_3. "adapt/fep"_fix_adapt_fep.html, "addtorque"_fix_addtorque.html, "atc"_fix_atc.html, "ave/correlate/long"_fix_ave_correlate_long.html, "colvars"_fix_colvars.html, "dpd/energy"_fix_dpd_energy.html, "drude"_fix_drude.html, "drude/transform/direct"_fix_drude_transform.html, "drude/transform/reverse"_fix_drude_transform.html, "eos/cv"_fix_eos_cv.html, "eos/table"_fix_eos_table.html, "eos/table/rx"_fix_eos_table_rx.html, "flow/gauss"_fix_flow_gauss.html, "gle"_fix_gle.html, "grem"_fix_grem.html, "imd"_fix_imd.html, "ipi"_fix_ipi.html, "langevin/drude"_fix_langevin_drude.html, "langevin/eff"_fix_langevin_eff.html, "lb/fluid"_fix_lb_fluid.html, "lb/momentum"_fix_lb_momentum.html, "lb/pc"_fix_lb_pc.html, "lb/rigid/pc/sphere"_fix_lb_rigid_pc_sphere.html, "lb/viscous"_fix_lb_viscous.html, "meso"_fix_meso.html, "manifoldforce"_fix_manifoldforce.html, "meso/stationary"_fix_meso_stationary.html, "nve/dot"_fix_nve_dot.html, "nve/dotc/langevin"_fix_nve_dotc_langevin.html, "nve/manifold/rattle"_fix_nve_manifold_rattle.html, "nvk"_fix_nvk.html, "nvt/manifold/rattle"_fix_nvt_manifold_rattle.html, "nph/eff"_fix_nh_eff.html, "npt/eff"_fix_nh_eff.html, "nve/eff"_fix_nve_eff.html, "nvt/eff"_fix_nh_eff.html, "nvt/sllod/eff"_fix_nvt_sllod_eff.html, "phonon"_fix_phonon.html, "pimd"_fix_pimd.html, "qbmsst"_fix_qbmsst.html, "qeq/reax"_fix_qeq_reax.html, "qmmm"_fix_qmmm.html, "qtb"_fix_qtb.html, "reax/c/bonds"_fix_reax_bonds.html, "reax/c/species"_fix_reaxc_species.html, "rx"_fix_rx.html, "saed/vtk"_fix_saed_vtk.html, "shardlow"_fix_shardlow.html, "smd"_fix_smd.html, "smd/adjust/dt"_fix_smd_adjust_dt.html, "smd/integrate/tlsph"_fix_smd_integrate_tlsph.html, "smd/integrate/ulsph"_fix_smd_integrate_ulsph.html, "smd/move/triangulated/surface"_fix_smd_move_triangulated_surface.html, "smd/setvel"_fix_smd_setvel.html, "smd/wall/surface"_fix_smd_wall_surface.html, "temp/rescale/eff"_fix_temp_rescale_eff.html, "ti/spring"_fix_ti_spring.html, "ttm/mod"_fix_ttm.html :tb(c=6,ea=c) :line Compute styles :h4 See the "compute"_compute.html command for one-line descriptions of each style or click on the style itself for a full description. Some of the styles have accelerated versions, which can be used if LAMMPS is built with the "appropriate accelerated package"_Section_accelerate.html. This is indicated by additional letters in parenthesis: g = GPU, i = USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT. "angle"_compute_angle.html, "angle/local"_compute_angle_local.html, "angmom/chunk"_compute_angmom_chunk.html, "body/local"_compute_body_local.html, "bond"_compute_bond.html, "bond/local"_compute_bond_local.html, "centro/atom"_compute_centro_atom.html, "chunk/atom"_compute_chunk_atom.html, "cluster/atom"_compute_cluster_atom.html, "cna/atom"_compute_cna_atom.html, "com"_compute_com.html, "com/chunk"_compute_com_chunk.html, "contact/atom"_compute_contact_atom.html, "coord/atom"_compute_coord_atom.html, "damage/atom"_compute_damage_atom.html, "dihedral"_compute_dihedral.html, "dihedral/local"_compute_dihedral_local.html, "dilatation/atom"_compute_dilatation_atom.html, "dipole/chunk"_compute_dipole_chunk.html, "displace/atom"_compute_displace_atom.html, "erotate/asphere"_compute_erotate_asphere.html, "erotate/rigid"_compute_erotate_rigid.html, "erotate/sphere"_compute_erotate_sphere.html, "erotate/sphere/atom"_compute_erotate_sphere_atom.html, "event/displace"_compute_event_displace.html, "global/atom"_compute_global_atom.html, "group/group"_compute_group_group.html, "gyration"_compute_gyration.html, "gyration/chunk"_compute_gyration_chunk.html, "heat/flux"_compute_heat_flux.html, "hexorder/atom"_compute_hexorder_atom.html, "improper"_compute_improper.html, "improper/local"_compute_improper_local.html, "inertia/chunk"_compute_inertia_chunk.html, "ke"_compute_ke.html, "ke/atom"_compute_ke_atom.html, "ke/rigid"_compute_ke_rigid.html, "msd"_compute_msd.html, "msd/chunk"_compute_msd_chunk.html, "msd/nongauss"_compute_msd_nongauss.html, "omega/chunk"_compute_omega_chunk.html, "orientorder/atom"_compute_orientorder_atom.html, "pair"_compute_pair.html, "pair/local"_compute_pair_local.html, "pe"_compute_pe.html, "pe/atom"_compute_pe_atom.html, "plasticity/atom"_compute_plasticity_atom.html, "pressure"_compute_pressure.html, "property/atom"_compute_property_atom.html, "property/local"_compute_property_local.html, "property/chunk"_compute_property_chunk.html, "rdf"_compute_rdf.html, "reduce"_compute_reduce.html, "reduce/region"_compute_reduce.html, "rigid/local"_compute_rigid_local.html, "slice"_compute_slice.html, "sna/atom"_compute_sna_atom.html, "snad/atom"_compute_sna_atom.html, "snav/atom"_compute_sna_atom.html, "stress/atom"_compute_stress_atom.html, "temp (k)"_compute_temp.html, "temp/asphere"_compute_temp_asphere.html, "temp/body"_compute_temp_body.html, "temp/chunk"_compute_temp_chunk.html, "temp/com"_compute_temp_com.html, "temp/deform"_compute_temp_deform.html, "temp/partial"_compute_temp_partial.html, "temp/profile"_compute_temp_profile.html, "temp/ramp"_compute_temp_ramp.html, "temp/region"_compute_temp_region.html, "temp/sphere"_compute_temp_sphere.html, "ti"_compute_ti.html, "torque/chunk"_compute_torque_chunk.html, "vacf"_compute_vacf.html, "vcm/chunk"_compute_vcm_chunk.html, "voronoi/atom"_compute_voronoi_atom.html :tb(c=6,ea=c) These are additional compute styles in USER packages, which can be used if "LAMMPS is built with the appropriate package"_Section_start.html#start_3. "ackland/atom"_compute_ackland_atom.html, "basal/atom"_compute_basal_atom.html, "dpd"_compute_dpd.html, "dpd/atom"_compute_dpd_atom.html, "fep"_compute_fep.html, "force/tally"_compute_tally.html, "heat/flux/tally"_compute_tally.html, "ke/eff"_compute_ke_eff.html, "ke/atom/eff"_compute_ke_atom_eff.html, "meso/e/atom"_compute_meso_e_atom.html, "meso/rho/atom"_compute_meso_rho_atom.html, "meso/t/atom"_compute_meso_t_atom.html, "pe/tally"_compute_tally.html, "pe/mol/tally"_compute_tally.html, "saed"_compute_saed.html, "smd/contact/radius"_compute_smd_contact_radius.html, "smd/damage"_compute_smd_damage.html, "smd/hourglass/error"_compute_smd_hourglass_error.html, "smd/internal/energy"_compute_smd_internal_energy.html, "smd/plastic/strain"_compute_smd_plastic_strain.html, "smd/plastic/strain/rate"_compute_smd_plastic_strain_rate.html, "smd/rho"_compute_smd_rho.html, "smd/tlsph/defgrad"_compute_smd_tlsph_defgrad.html, "smd/tlsph/dt"_compute_smd_tlsph_dt.html, "smd/tlsph/num/neighs"_compute_smd_tlsph_num_neighs.html, "smd/tlsph/shape"_compute_smd_tlsph_shape.html, "smd/tlsph/strain"_compute_smd_tlsph_strain.html, "smd/tlsph/strain/rate"_compute_smd_tlsph_strain_rate.html, "smd/tlsph/stress"_compute_smd_tlsph_stress.html, "smd/triangle/mesh/vertices"_compute_smd_triangle_mesh_vertices.html, "smd/ulsph/num/neighs"_compute_smd_ulsph_num_neighs.html, "smd/ulsph/strain"_compute_smd_ulsph_strain.html, "smd/ulsph/strain/rate"_compute_smd_ulsph_strain_rate.html, "smd/ulsph/stress"_compute_smd_ulsph_stress.html, "smd/vol"_compute_smd_vol.html, "stress/tally"_compute_tally.html, "temp/drude"_compute_temp_drude.html, "temp/eff"_compute_temp_eff.html, "temp/deform/eff"_compute_temp_deform_eff.html, "temp/region/eff"_compute_temp_region_eff.html, "temp/rotate"_compute_temp_rotate.html, "xrd"_compute_xrd.html :tb(c=6,ea=c) :line Pair_style potentials :h4 See the "pair_style"_pair_style.html command for an overview of pair potentials. Click on the style itself for a full description. Many of the styles have accelerated versions, which can be used if LAMMPS is built with the "appropriate accelerated package"_Section_accelerate.html. This is indicated by additional letters in parenthesis: g = GPU, i = USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT. "none"_pair_none.html, "zero"_pair_zero.html, "hybrid"_pair_hybrid.html, "hybrid/overlay"_pair_hybrid.html, "adp (o)"_pair_adp.html, "airebo (o)"_pair_airebo.html, "airebo/morse (o)"_pair_airebo.html, "beck (go)"_pair_beck.html, "body"_pair_body.html, "bop"_pair_bop.html, "born (go)"_pair_born.html, "born/coul/dsf"_pair_born.html, "born/coul/dsf/cs"_pair_born.html, "born/coul/long (go)"_pair_born.html, "born/coul/long/cs"_pair_born.html, "born/coul/msm (o)"_pair_born.html, "born/coul/wolf (go)"_pair_born.html, "brownian (o)"_pair_brownian.html, "brownian/poly (o)"_pair_brownian.html, "buck (gkio)"_pair_buck.html, "buck/coul/cut (gkio)"_pair_buck.html, "buck/coul/long (gkio)"_pair_buck.html, "buck/coul/long/cs"_pair_buck.html, "buck/coul/msm (o)"_pair_buck.html, "buck/long/coul/long (o)"_pair_buck_long.html, "colloid (go)"_pair_colloid.html, "comb (o)"_pair_comb.html, "comb3"_pair_comb.html, "coul/cut (gko)"_pair_coul.html, "coul/debye (gko)"_pair_coul.html, "coul/dsf (gko)"_pair_coul.html, "coul/long (gko)"_pair_coul.html, "coul/long/cs"_pair_coul.html, "coul/msm"_pair_coul.html, "coul/streitz"_pair_coul.html, "coul/wolf (ko)"_pair_coul.html, "dpd (go)"_pair_dpd.html, "dpd/tstat (go)"_pair_dpd.html, "dsmc"_pair_dsmc.html, "eam (gkiot)"_pair_eam.html, "eam/alloy (gkot)"_pair_eam.html, "eam/fs (gkot)"_pair_eam.html, "eim (o)"_pair_eim.html, "gauss (go)"_pair_gauss.html, "gayberne (gio)"_pair_gayberne.html, "gran/hertz/history (o)"_pair_gran.html, "gran/hooke (o)"_pair_gran.html, "gran/hooke/history (o)"_pair_gran.html, "hbond/dreiding/lj (o)"_pair_hbond_dreiding.html, "hbond/dreiding/morse (o)"_pair_hbond_dreiding.html, "kim"_pair_kim.html, "lcbop"_pair_lcbop.html, "line/lj"_pair_line_lj.html, "lj/charmm/coul/charmm (ko)"_pair_charmm.html, "lj/charmm/coul/charmm/implicit (ko)"_pair_charmm.html, "lj/charmm/coul/long (giko)"_pair_charmm.html, "lj/charmm/coul/msm"_pair_charmm.html, "lj/class2 (gko)"_pair_class2.html, "lj/class2/coul/cut (ko)"_pair_class2.html, "lj/class2/coul/long (gko)"_pair_class2.html, "lj/cubic (go)"_pair_lj_cubic.html, "lj/cut (gikot)"_pair_lj.html, "lj/cut/coul/cut (gko)"_pair_lj.html, "lj/cut/coul/debye (gko)"_pair_lj.html, "lj/cut/coul/dsf (gko)"_pair_lj.html, "lj/cut/coul/long (gikot)"_pair_lj.html, "lj/cut/coul/long/cs"_pair_lj.html, "lj/cut/coul/msm (go)"_pair_lj.html, "lj/cut/dipole/cut (go)"_pair_dipole.html, "lj/cut/dipole/long"_pair_dipole.html, "lj/cut/tip4p/cut (o)"_pair_lj.html, "lj/cut/tip4p/long (ot)"_pair_lj.html, "lj/expand (gko)"_pair_lj_expand.html, "lj/gromacs (gko)"_pair_gromacs.html, "lj/gromacs/coul/gromacs (ko)"_pair_gromacs.html, "lj/long/coul/long (o)"_pair_lj_long.html, "lj/long/dipole/long"_pair_dipole.html, "lj/long/tip4p/long"_pair_lj_long.html, "lj/smooth (o)"_pair_lj_smooth.html, "lj/smooth/linear (o)"_pair_lj_smooth_linear.html, "lj96/cut (go)"_pair_lj96.html, "lubricate (o)"_pair_lubricate.html, "lubricate/poly (o)"_pair_lubricate.html, "lubricateU"_pair_lubricateU.html, "lubricateU/poly"_pair_lubricateU.html, "meam"_pair_meam.html, "mie/cut (o)"_pair_mie.html, -"morse (got)"_pair_morse.html, +"morse (gkot)"_pair_morse.html, "nb3b/harmonic (o)"_pair_nb3b_harmonic.html, "nm/cut (o)"_pair_nm.html, "nm/cut/coul/cut (o)"_pair_nm.html, "nm/cut/coul/long (o)"_pair_nm.html, "peri/eps"_pair_peri.html, "peri/lps (o)"_pair_peri.html, "peri/pmb (o)"_pair_peri.html, "peri/ves"_pair_peri.html, "polymorphic"_pair_polymorphic.html, "reax"_pair_reax.html, "rebo (o)"_pair_airebo.html, "resquared (go)"_pair_resquared.html, "snap"_pair_snap.html, "soft (go)"_pair_soft.html, "sw (gkio)"_pair_sw.html, "table (gko)"_pair_table.html, "tersoff (gkio)"_pair_tersoff.html, "tersoff/mod (gko)"_pair_tersoff_mod.html, "tersoff/mod/c (o)"_pair_tersoff_mod.html, "tersoff/zbl (gko)"_pair_tersoff_zbl.html, "tip4p/cut (o)"_pair_coul.html, "tip4p/long (o)"_pair_coul.html, "tri/lj"_pair_tri_lj.html, "vashishta (ko)"_pair_vashishta.html, "vashishta/table (o)"_pair_vashishta.html, "yukawa (go)"_pair_yukawa.html, "yukawa/colloid (go)"_pair_yukawa_colloid.html, "zbl (go)"_pair_zbl.html :tb(c=4,ea=c) These are additional pair styles in USER packages, which can be used if "LAMMPS is built with the appropriate package"_Section_start.html#start_3. "agni (o)"_pair_agni.html, "awpmd/cut"_pair_awpmd.html, "buck/mdf"_pair_mdf.html, "coul/cut/soft (o)"_pair_lj_soft.html, "coul/diel (o)"_pair_coul_diel.html, "coul/long/soft (o)"_pair_lj_soft.html, "dpd/fdt"_pair_dpd_fdt.html, "dpd/fdt/energy"_pair_dpd_fdt.html, "eam/cd (o)"_pair_eam.html, "edip (o)"_pair_edip.html, "eff/cut"_pair_eff.html, "exp6/rx"_pair_exp6_rx.html, "gauss/cut"_pair_gauss.html, "lennard/mdf"_pair_mdf.html, "list"_pair_list.html, "lj/charmm/coul/long/soft (o)"_pair_charmm.html, "lj/cut/coul/cut/soft (o)"_pair_lj_soft.html, "lj/cut/coul/long/soft (o)"_pair_lj_soft.html, "lj/cut/dipole/sf (go)"_pair_dipole.html, "lj/cut/soft (o)"_pair_lj_soft.html, "lj/cut/thole/long (o)"_pair_thole.html, "lj/cut/tip4p/long/soft (o)"_pair_lj_soft.html, "lj/mdf"_pair_mdf.html, "lj/sdk (gko)"_pair_sdk.html, "lj/sdk/coul/long (go)"_pair_sdk.html, "lj/sdk/coul/msm (o)"_pair_sdk.html, "lj/sf (o)"_pair_lj_sf.html, "meam/spline (o)"_pair_meam_spline.html, "meam/sw/spline"_pair_meam_sw_spline.html, "mgpt"_pair_mgpt.html, "morse/smooth/linear"_pair_morse.html, "morse/soft"_pair_morse.html, "multi/lucy"_pair_multi_lucy.html, "multi/lucy/rx"_pair_multi_lucy_rx.html, "oxdna/coaxstk"_pair_oxdna.html, "oxdna/excv"_pair_oxdna.html, "oxdna/hbond"_pair_oxdna.html, "oxdna/stk"_pair_oxdna.html, "oxdna/xstk"_pair_oxdna.html, "quip"_pair_quip.html, "reax/c (k)"_pair_reax_c.html, "smd/hertz"_pair_smd_hertz.html, "smd/tlsph"_pair_smd_tlsph.html, "smd/triangulated/surface"_pair_smd_triangulated_surface.html, "smd/ulsph"_pair_smd_ulsph.html, "smtbq"_pair_smtbq.html, "sph/heatconduction"_pair_sph_heatconduction.html, "sph/idealgas"_pair_sph_idealgas.html, "sph/lj"_pair_sph_lj.html, "sph/rhosum"_pair_sph_rhosum.html, "sph/taitwater"_pair_sph_taitwater.html, "sph/taitwater/morris"_pair_sph_taitwater_morris.html, "srp"_pair_srp.html, "table/rx"_pair_table_rx.html, "tersoff/table (o)"_pair_tersoff.html, "thole"_pair_thole.html, "tip4p/long/soft (o)"_pair_lj_soft.html :tb(c=4,ea=c) :line Bond_style potentials :h4 See the "bond_style"_bond_style.html command for an overview of bond potentials. Click on the style itself for a full description. Some of the styles have accelerated versions, which can be used if LAMMPS is built with the "appropriate accelerated package"_Section_accelerate.html. This is indicated by additional letters in parenthesis: g = GPU, i = USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT. "none"_bond_none.html, "zero"_bond_zero.html, "hybrid"_bond_hybrid.html, "class2 (o)"_bond_class2.html, "fene (iko)"_bond_fene.html, "fene/expand (o)"_bond_fene_expand.html, "harmonic (ko)"_bond_harmonic.html, "morse (o)"_bond_morse.html, "nonlinear (o)"_bond_nonlinear.html, "quartic (o)"_bond_quartic.html, "table (o)"_bond_table.html :tb(c=4,ea=c) These are additional bond styles in USER packages, which can be used if "LAMMPS is built with the appropriate package"_Section_start.html#start_3. "harmonic/shift (o)"_bond_harmonic_shift.html, "harmonic/shift/cut (o)"_bond_harmonic_shift_cut.html, "oxdna/fene"_bond_oxdna_fene.html :tb(c=4,ea=c) :line Angle_style potentials :h4 See the "angle_style"_angle_style.html command for an overview of angle potentials. Click on the style itself for a full description. Some of the styles have accelerated versions, which can be used if LAMMPS is built with the "appropriate accelerated package"_Section_accelerate.html. This is indicated by additional letters in parenthesis: g = GPU, i = USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT. "none"_angle_none.html, "zero"_angle_zero.html, "hybrid"_angle_hybrid.html, "charmm (ko)"_angle_charmm.html, "class2 (o)"_angle_class2.html, "cosine (o)"_angle_cosine.html, "cosine/delta (o)"_angle_cosine_delta.html, "cosine/periodic (o)"_angle_cosine_periodic.html, "cosine/squared (o)"_angle_cosine_squared.html, "harmonic (iko)"_angle_harmonic.html, "table (o)"_angle_table.html :tb(c=4,ea=c) These are additional angle styles in USER packages, which can be used if "LAMMPS is built with the appropriate package"_Section_start.html#start_3. "cosine/shift (o)"_angle_cosine_shift.html, "cosine/shift/exp (o)"_angle_cosine_shift_exp.html, "dipole (o)"_angle_dipole.html, "fourier (o)"_angle_fourier.html, "fourier/simple (o)"_angle_fourier_simple.html, "quartic (o)"_angle_quartic.html, "sdk"_angle_sdk.html :tb(c=4,ea=c) :line Dihedral_style potentials :h4 See the "dihedral_style"_dihedral_style.html command for an overview of dihedral potentials. Click on the style itself for a full description. Some of the styles have accelerated versions, which can be used if LAMMPS is built with the "appropriate accelerated package"_Section_accelerate.html. This is indicated by additional letters in parenthesis: g = GPU, i = USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT. "none"_dihedral_none.html, "zero"_dihedral_zero.html, "hybrid"_dihedral_hybrid.html, "charmm (ko)"_dihedral_charmm.html, "class2 (o)"_dihedral_class2.html, "harmonic (io)"_dihedral_harmonic.html, "helix (o)"_dihedral_helix.html, "multi/harmonic (o)"_dihedral_multi_harmonic.html, "opls (iko)"_dihedral_opls.html :tb(c=4,ea=c) These are additional dihedral styles in USER packages, which can be used if "LAMMPS is built with the appropriate package"_Section_start.html#start_3. "cosine/shift/exp (o)"_dihedral_cosine_shift_exp.html, "fourier (o)"_dihedral_fourier.html, "nharmonic (o)"_dihedral_nharmonic.html, "quadratic (o)"_dihedral_quadratic.html, "spherical (o)"_dihedral_spherical.html, "table (o)"_dihedral_table.html :tb(c=4,ea=c) :line Improper_style potentials :h4 See the "improper_style"_improper_style.html command for an overview of improper potentials. Click on the style itself for a full description. Some of the styles have accelerated versions, which can be used if LAMMPS is built with the "appropriate accelerated package"_Section_accelerate.html. This is indicated by additional letters in parenthesis: g = GPU, i = USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT. "none"_improper_none.html, "zero"_improper_zero.html, "hybrid"_improper_hybrid.html, "class2 (o)"_improper_class2.html, "cvff (io)"_improper_cvff.html, "harmonic (ko)"_improper_harmonic.html, "umbrella (o)"_improper_umbrella.html :tb(c=4,ea=c) These are additional improper styles in USER packages, which can be used if "LAMMPS is built with the appropriate package"_Section_start.html#start_3. "cossq (o)"_improper_cossq.html, "distance"_improper_distance.html, "fourier (o)"_improper_fourier.html, "ring (o)"_improper_ring.html :tb(c=4,ea=c) :line Kspace solvers :h4 See the "kspace_style"_kspace_style.html command for an overview of Kspace solvers. Click on the style itself for a full description. Some of the styles have accelerated versions, which can be used if LAMMPS is built with the "appropriate accelerated package"_Section_accelerate.html. This is indicated by additional letters in parenthesis: g = GPU, i = USER-INTEL, k = KOKKOS, o = USER-OMP, t = OPT. "ewald (o)"_kspace_style.html, "ewald/disp"_kspace_style.html, "msm (o)"_kspace_style.html, "msm/cg (o)"_kspace_style.html, "pppm (go)"_kspace_style.html, "pppm/cg (o)"_kspace_style.html, "pppm/disp"_kspace_style.html, "pppm/disp/tip4p"_kspace_style.html, "pppm/stagger"_kspace_style.html, "pppm/tip4p (o)"_kspace_style.html :tb(c=4,ea=c) diff --git a/doc/src/compute_group_group.txt b/doc/src/compute_group_group.txt index 3ea229cca..175ada957 100644 --- a/doc/src/compute_group_group.txt +++ b/doc/src/compute_group_group.txt @@ -1,127 +1,156 @@ "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 group/group command :h3 [Syntax:] compute ID group-ID group/group group2-ID keyword value ... :pre ID, group-ID are documented in "compute"_compute.html command :ulb,l group/group = style name of this compute command :l group2-ID = group ID of second (or same) group :l zero or more keyword/value pairs may be appended :l keyword = {pair} or {kspace} or {boundary} :l {pair} value = {yes} or {no} {kspace} value = {yes} or {no} {boundary} value = {yes} or {no} :pre + {molecule} value = {off} or {inter} or {intra} :pre :ule [Examples:] compute 1 lower group/group upper compute 1 lower group/group upper kspace yes compute mine fluid group/group wall :pre [Description:] Define a computation that calculates the total energy and force interaction between two groups of atoms: the compute group and the specified group2. The two groups can be the same. If the {pair} keyword is set to {yes}, which is the default, then the the interaction energy will include a pair component which is defined as the pairwise energy between all pairs of atoms where one atom in the pair is in the first group and the other is in the second group. Likewise, the interaction force calculated by this compute will include the force on the compute group atoms due to pairwise interactions with atoms in the specified group2. NOTE: The energies computed by the {pair} keyword do not include tail corrections, even if they are enabled via the "pair_modify"_pair_modify.html command. +If the {molecule} keyword is set to {inter} or {intra} than an +additional check is made based on the molecule IDs of the two atoms in +each pair before including their pairwise interaction energy and +force. For the {inter} setting, the two atoms must be in different +molecules. For the {intra} setting, the two atoms must be in the same +molecule. + If the {kspace} keyword is set to {yes}, which is not the default, and if a "kspace_style"_kspace_style.html is defined, then the interaction energy will include a Kspace component which is the long-range Coulombic energy between all the atoms in the first group and all the atoms in the 2nd group. Likewise, the interaction force calculated by this compute will include the force on the compute group atoms due to long-range Coulombic interactions with atoms in the specified group2. Normally the long-range Coulombic energy converges only when the net charge of the unit cell is zero. However, one can assume the net charge of the system is neutralized by a uniform background plasma, and a correction to the system energy can be applied to reduce artifacts. For more information see "(Bogusz)"_#Bogusz. If the {boundary} keyword is set to {yes}, which is the default, and {kspace} contributions are included, then this energy correction term will be added to the total group-group energy. This correction term does not affect the force calculation and will be zero if one or both of the groups are charge neutral. This energy correction term is the same as that included in the regular Ewald and PPPM routines. +NOTE: The {molecule} setting only affects the group/group +contributions calculated by the {pair} keyword. It does not affect +the group/group contributions calculated by the {kspace} keyword. + This compute does not calculate any bond or angle or dihedral or improper interactions between atoms in the two groups. :line The pairwise contributions to the group-group interactions are calculated by looping over a neighbor list. The Kspace contribution to the group-group interactions require essentially the same amount of work (FFTs, Ewald summation) as computing long-range forces for the entire system. Thus it can be costly to invoke this compute too frequently. +NOTE: If you have a bonded system, then the settings of +"special_bonds"_special_bonds.html command can remove pairwise +interactions between atoms in the same bond, angle, or dihedral. This +is the default setting for the "special_bonds"_special_bonds.html +command, and means those pairwise interactions do not appear in the +neighbor list. Because this compute uses a neighbor list, it also +means those pairs will not be included in the group/group interaction. +This does not apply when using long-range coulomb interactions +({coul/long}, {coul/msm}, {coul/wolf} or similar. One way to get +around this would be to set special_bond scaling factors to very tiny +numbers that are not exactly zero (e.g. 1.0e-50). Another workaround +is to write a dump file, and use the "rerun"_rerun.html command to +compute the group/group interactions for snapshots in the dump file. +The rerun script can use a "special_bonds"_special_bonds.html command +that includes all pairs in the neighbor list. + If you desire a breakdown of the interactions into a pairwise and Kspace component, simply invoke the compute twice with the appropriate yes/no settings for the {pair} and {kspace} keywords. This is no more costly than using a single compute with both keywords set to {yes}. The individual contributions can be summed in a "variable"_variable.html if desired. This "document"_PDF/kspace.pdf describes how the long-range group-group calculations are performed. :line [Output info:] This compute calculates a global scalar (the energy) and a global vector of length 3 (force), which can be accessed by indices 1-3. These values can be used by any command that uses global scalar or vector values from a compute as input. See "this section"_Section_howto.html#howto_15 for an overview of LAMMPS output options. Both the scalar and vector values calculated by this compute are "extensive". The scalar value will be in energy "units"_units.html. The vector values will be in force "units"_units.html. [Restrictions:] Not all pair styles can be evaluated in a pairwise mode as required by this compute. For example, 3-body and other many-body potentials, such as "Tersoff"_pair_tersoff.html and "Stillinger-Weber"_pair_sw.html cannot be used. "EAM"_pair_eam.html potentials only include the pair potential portion of the EAM interaction when used by this compute, not the embedding term. Not all Kspace styles support calculation of group/group interactions. The {ewald} and {pppm} styles do. [Related commands:] none [Default:] -The option defaults are pair = yes, kspace = no, and boundary = yes. +The option defaults are pair = yes, kspace = no, boundary = yes, +molecule = off. :line :link(Bogusz) Bogusz et al, J Chem Phys, 108, 7070 (1998) diff --git a/examples/VISCOSITY/in.einstein.2d b/examples/VISCOSITY/in.einstein.2d index 72ccf5d7a..d728ee620 100644 --- a/examples/VISCOSITY/in.einstein.2d +++ b/examples/VISCOSITY/in.einstein.2d @@ -1,84 +1,85 @@ # sample LAMMPS input script for viscosity of 2d LJ liquid # Einstein form of Green-Kubo # settings variable x equal 20 variable y equal 20 variable rho equal 0.6 variable t equal 1.0 variable rc equal 2.5 variable p equal 400 # correlation length variable s equal 5 # sample interval variable d equal $p*$s # dump interval # problem setup units lj dimension 2 atom_style atomic neigh_modify delay 0 every 1 lattice sq2 ${rho} region simbox block 0 $x 0 $y -0.1 0.1 create_box 1 simbox create_atoms 1 box pair_style lj/cut ${rc} pair_coeff * * 1 1 mass * 1.0 velocity all create $t 97287 # equilibration run fix 1 all nve fix 2 all langevin $t $t 0.1 498094 fix 3 all enforce2d thermo $d run 10000 velocity all scale $t unfix 2 # Einstein viscosity calculation reset_timestep 0 # Define distinct components of symmetric traceless stress tensor variable pxy equal pxy variable pxx equal pxx-press -fix avstress all ave/time $s $p $d v_pxy v_pxx ave one file einstein.dat +fix avstress all ave/time $s $p $d v_pxy v_pxx ave one & + file profile.einstein.2d # Diagonal components of SS are larger by factor 2-2/d, # which is 4/3 for d=3, but 1 for d=2. # See Daivis and Evans, J.Chem.Phys, 100, 541-547 (1994) variable scale equal vol/(2.0*$t*dt*$d) variable diagfac equal 2-2/2 variable deltasqxy equal (f_avstress[1]*$d*dt)^2 variable deltasqxx equal (f_avstress[2]*$d*dt)^2/${diagfac} # compute mean square displacements as running averages fix avdeltasq all ave/time $d 1 $d v_deltasqxy v_deltasqxx ave running # convert to viscosities variable vxy equal f_avdeltasq[1]*${scale} variable vxx equal f_avdeltasq[2]*${scale} thermo_style custom step temp pe press pxy v_vxy v_vxx run 500000 variable etaxy equal v_vxy variable etaxx equal v_vxx variable eta equal 0.5*(${etaxy}+${etaxx}) print "running average viscosity: ${eta}" diff --git a/examples/coreshell/log.9Nov16.coreshell.dsf.1 b/examples/coreshell/log.9Nov16.coreshell.dsf.g++.1 similarity index 100% rename from examples/coreshell/log.9Nov16.coreshell.dsf.1 rename to examples/coreshell/log.9Nov16.coreshell.dsf.g++.1 diff --git a/examples/coreshell/log.9Nov16.coreshell.dsf.g++.4 b/examples/coreshell/log.9Nov16.coreshell.dsf.g++.4 new file mode 100644 index 000000000..ccf1b4b03 --- /dev/null +++ b/examples/coreshell/log.9Nov16.coreshell.dsf.g++.4 @@ -0,0 +1,189 @@ +LAMMPS (26 Jan 2017) +# Testsystem for core-shell model compared to Mitchel and Finchham +# Hendrik Heenen, June 2014 + +# ------------------------ INITIALIZATION ---------------------------- + +units metal +dimension 3 +boundary p p p +atom_style full + +# ----------------------- ATOM DEFINITION ---------------------------- + +fix csinfo all property/atom i_CSID +read_data data.coreshell fix csinfo NULL CS-Info + orthogonal box = (0 0 0) to (24.096 24.096 24.096) + 1 by 2 by 2 MPI processor grid + reading atoms ... + 432 atoms + scanning bonds ... + 1 = max bonds/atom + reading bonds ... + 216 bonds + 1 = max # of 1-2 neighbors + 0 = max # of 1-3 neighbors + 0 = max # of 1-4 neighbors + 1 = max # of special neighbors + +group cores type 1 2 +216 atoms in group cores +group shells type 3 4 +216 atoms in group shells + +neighbor 2.0 bin +comm_modify vel yes + +# ------------------------ FORCE FIELDS ------------------------------ + +pair_style born/coul/dsf/cs 0.1 20.0 20.0 # A, rho, sigma=0, C, D +pair_coeff * * 0.0 1.000 0.00 0.00 0.00 +pair_coeff 3 3 487.0 0.23768 0.00 1.05 0.50 #Na-Na +pair_coeff 3 4 145134.0 0.23768 0.00 6.99 8.70 #Na-Cl +pair_coeff 4 4 405774.0 0.23768 0.00 72.40 145.40 #Cl-Cl + +bond_style harmonic +bond_coeff 1 63.014 0.0 +bond_coeff 2 25.724 0.0 + +# ------------------------ Equilibration Run ------------------------------- + +reset_timestep 0 + +thermo 50 +thermo_style custom step etotal pe ke temp press epair evdwl ecoul elong ebond fnorm fmax vol + +compute CSequ all temp/cs cores shells + +# output via chunk method + +#compute prop all property/atom i_CSID +#compute cs_chunk all chunk/atom c_prop +#compute cstherm all temp/chunk cs_chunk temp internal com yes cdof 3.0 +#fix ave_chunk all ave/time 100 1 100 c_cstherm file chunk.dump mode vector + +thermo_modify temp CSequ + +# velocity bias option + +velocity all create 1427 134 dist gaussian mom yes rot no bias yes temp CSequ +Neighbor list info ... + update every 1 steps, delay 10 steps, check yes + max neighbors/atom: 2000, page size: 100000 + master list distance cutoff = 22 + ghost atom cutoff = 22 + binsize = 11, bins = 3 3 3 + 1 neighbor lists, perpetual/occasional/extra = 1 0 0 + (1) pair born/coul/dsf/cs, half, perpetual + pair build: half/bin/newton + stencil: half/bin/3d/newton + bin: standard +velocity all scale 1427 temp CSequ + +fix thermoberendsen all temp/berendsen 1427 1427 0.4 +fix nve all nve +fix_modify thermoberendsen temp CSequ + +# 2 fmsec timestep + +timestep 0.002 +run 500 +Memory usage per processor = 6.8559 Mbytes +Step TotEng PotEng KinEng Temp Press E_pair E_vdwl E_coul E_long E_bond Fnorm Fmax Volume + 0 -635.80596 -675.46362 39.657659 1427 -21302.622 -675.46362 1.6320365 -677.09565 0 0 1.5814015e-14 3.2317898e-15 13990.5 + 50 -634.07021 -666.11867 32.048452 1153.1982 -4560.945 -668.28236 37.756542 -706.0389 0 2.163691 13.802484 3.022372 13990.5 + 100 -631.97128 -662.02544 30.054164 1081.4378 -3497.564 -664.61825 39.275003 -703.89325 0 2.5928078 13.956833 2.5417699 13990.5 + 150 -630.14953 -663.04215 32.892622 1183.5739 -88.43828 -665.63444 46.239965 -711.87441 0 2.5922927 14.667898 2.4964255 13990.5 + 200 -628.52878 -663.9795 35.45072 1275.6219 -1755.9004 -666.73564 41.758052 -708.49369 0 2.7561421 14.230743 3.0924004 13990.5 + 250 -627.27102 -662.025 34.753978 1250.5511 -1234.0918 -665.13519 43.170874 -708.30606 0 3.1101887 14.221086 1.941354 13990.5 + 300 -626.5495 -663.74287 37.193368 1338.3275 -2049.3444 -666.45574 40.476148 -706.93188 0 2.7128711 13.330425 1.7756755 13990.5 + 350 -625.87313 -665.21855 39.345421 1415.7647 -1543.1723 -667.90872 41.577366 -709.48609 0 2.6901682 13.541311 1.854662 13990.5 + 400 -625.09344 -661.26404 36.1706 1301.5253 -729.96729 -664.10334 43.468765 -707.57211 0 2.8392963 13.663555 1.9067551 13990.5 + 450 -624.46214 -660.01362 35.551477 1279.2474 -1617.7158 -663.06571 41.644856 -704.71057 0 3.0520921 14.527005 1.7280213 13990.5 + 500 -623.49246 -659.2527 35.76024 1286.7593 -935.99238 -662.32953 43.038808 -705.36834 0 3.0768302 14.099593 1.9831106 13990.5 +Loop time of 4.09864 on 4 procs for 500 steps with 432 atoms + +Performance: 21.080 ns/day, 1.139 hours/ns, 121.992 timesteps/s +99.7% CPU use with 4 MPI tasks x no OpenMP threads + +MPI task timing breakdown: +Section | min time | avg time | max time |%varavg| %total +--------------------------------------------------------------- +Pair | 3.3804 | 3.568 | 3.8354 | 8.9 | 87.05 +Bond | 0.00074339 | 0.00079519 | 0.00087976 | 0.0 | 0.02 +Neigh | 0.045851 | 0.046084 | 0.046361 | 0.1 | 1.12 +Comm | 0.20413 | 0.47123 | 0.65875 | 24.3 | 11.50 +Output | 0.00044298 | 0.00046057 | 0.00051165 | 0.0 | 0.01 +Modify | 0.0064909 | 0.0067219 | 0.0069766 | 0.2 | 0.16 +Other | | 0.005345 | | | 0.13 + +Nlocal: 108 ave 114 max 105 min +Histogram: 1 1 1 0 0 0 0 0 0 1 +Nghost: 6527 ave 6599 max 6472 min +Histogram: 1 0 1 0 1 0 0 0 0 1 +Neighs: 74388.2 ave 75855 max 73680 min +Histogram: 1 2 0 0 0 0 0 0 0 1 + +Total # of neighbors = 297553 +Ave neighs/atom = 688.78 +Ave special neighs/atom = 1 +Neighbor list builds = 20 +Dangerous builds = 0 + +unfix thermoberendsen + +# ------------------------ Dynamic Run ------------------------------- + +run 1000 +Memory usage per processor = 6.85787 Mbytes +Step TotEng PotEng KinEng Temp Press E_pair E_vdwl E_coul E_long E_bond Fnorm Fmax Volume + 500 -623.49319 -659.2527 35.759511 1286.7331 -936.04802 -662.32953 43.038808 -705.36834 0 3.0768302 14.099593 1.9831106 13990.5 + 550 -623.44059 -663.57938 40.138795 1444.3127 -935.73484 -666.2789 42.563337 -708.84224 0 2.6995167 13.918509 2.3189805 13990.5 + 600 -623.4703 -660.01592 36.545618 1315.0196 1327.3492 -663.08845 47.985462 -711.07391 0 3.0725254 15.192713 2.4098428 13990.5 + 650 -623.46796 -661.56776 38.099807 1370.9439 457.82439 -664.81976 45.495622 -710.31538 0 3.2519966 15.026057 1.8500226 13990.5 + 700 -623.50158 -659.5131 36.011523 1295.8012 -460.03772 -663.1078 43.938203 -707.046 0 3.5946908 14.660979 2.4825518 13990.5 + 750 -623.44787 -661.93353 38.485658 1384.8279 97.429626 -664.9551 45.083146 -710.03825 0 3.0215753 15.10043 2.3433897 13990.5 + 800 -623.48215 -659.50655 36.024402 1296.2647 1097.3866 -662.61124 47.251998 -709.86324 0 3.1046914 14.556382 2.0543766 13990.5 + 850 -623.45868 -661.13782 37.679134 1355.8068 -1802.1624 -664.41257 40.70845 -705.12102 0 3.2747525 14.691444 2.2054332 13990.5 + 900 -623.43556 -663.59137 40.155815 1444.9251 534.99197 -666.71877 45.601619 -712.32039 0 3.127395 14.741411 2.5807895 13990.5 + 950 -623.51318 -661.57916 38.06598 1369.7267 -678.12625 -664.37535 43.207862 -707.58322 0 2.7961988 14.430307 2.3936105 13990.5 + 1000 -623.47287 -661.22274 37.749874 1358.3523 634.7979 -664.42973 46.373361 -710.80309 0 3.2069879 15.891192 2.4042765 13990.5 + 1050 -623.48133 -661.52868 38.047347 1369.0562 -583.15228 -664.6098 43.618772 -708.22857 0 3.081116 14.806856 2.3447613 13990.5 + 1100 -623.47867 -661.83761 38.358946 1380.2685 -868.9779 -664.8826 42.84846 -707.73106 0 3.044983 14.69567 2.399143 13990.5 + 1150 -623.44713 -661.21299 37.765857 1358.9274 405.14554 -664.09567 45.578739 -709.6744 0 2.8826753 15.437367 3.1381305 13990.5 + 1200 -623.46549 -660.91706 37.451568 1347.6183 699.78996 -664.0883 46.36297 -710.45127 0 3.1712473 15.109665 1.8891886 13990.5 + 1250 -623.49296 -658.2218 34.728838 1249.6464 1061.0154 -661.29052 47.668699 -708.95922 0 3.0687228 14.901367 2.3964137 13990.5 + 1300 -623.49837 -660.91022 37.411844 1346.1889 226.99512 -664.35989 45.352287 -709.71217 0 3.4496704 15.161542 2.2137993 13990.5 + 1350 -623.46718 -658.80365 35.336469 1271.5108 1039.6469 -662.16908 47.565671 -709.73475 0 3.3654314 15.892516 2.7888426 13990.5 + 1400 -623.47124 -661.45375 37.982513 1366.7233 -379.56023 -664.6321 43.788306 -708.42041 0 3.1783497 14.251126 1.7415409 13990.5 + 1450 -623.46671 -660.17518 36.708464 1320.8792 -374.37056 -662.92706 44.083648 -707.01071 0 2.7518803 15.210167 1.9984277 13990.5 + 1500 -623.50515 -659.06488 35.559725 1279.5442 260.37822 -662.39548 45.779764 -708.17524 0 3.3306005 14.682396 2.4201107 13990.5 +Loop time of 8.26746 on 4 procs for 1000 steps with 432 atoms + +Performance: 20.901 ns/day, 1.148 hours/ns, 120.956 timesteps/s +99.7% CPU use with 4 MPI tasks x no OpenMP threads + +MPI task timing breakdown: +Section | min time | avg time | max time |%varavg| %total +--------------------------------------------------------------- +Pair | 6.706 | 7.1568 | 7.6597 | 12.7 | 86.57 +Bond | 0.0014617 | 0.0015531 | 0.0016506 | 0.2 | 0.02 +Neigh | 0.10511 | 0.10522 | 0.10532 | 0.0 | 1.27 +Comm | 0.48547 | 0.98841 | 1.4393 | 34.0 | 11.96 +Output | 0.0012085 | 0.0012462 | 0.0013196 | 0.1 | 0.02 +Modify | 0.0021446 | 0.0021989 | 0.0022545 | 0.1 | 0.03 +Other | | 0.01204 | | | 0.15 + +Nlocal: 108 ave 114 max 94 min +Histogram: 1 0 0 0 0 0 0 0 1 2 +Nghost: 6512.25 ave 6586 max 6456 min +Histogram: 1 0 0 2 0 0 0 0 0 1 +Neighs: 74248.2 ave 77441 max 65858 min +Histogram: 1 0 0 0 0 0 0 0 0 3 + +Total # of neighbors = 296993 +Ave neighs/atom = 687.484 +Ave special neighs/atom = 1 +Neighbor list builds = 46 +Dangerous builds = 0 +Total wall time: 0:00:12 diff --git a/potentials/Ge.tersoff b/potentials/Ge.tersoff new file mode 100644 index 000000000..2904f3a63 --- /dev/null +++ b/potentials/Ge.tersoff @@ -0,0 +1,17 @@ +# DATE: 2016-12-21 CONTRIBUTOR: Sayyed Jalil Mahdizadh, saja.mahdizadeh@gmail.com CITATION: Mahdizadeh, Akhlamadi, J. Mol. Graph. Model. 72, 1-5 (2017) + +# Tersoff parameters for various elements and mixtures +# multiple entries can be added to this file, LAMMPS reads the ones it needs +# these entries are in LAMMPS "metal" units: +# A,B = eV; lambda1,lambda2,lambda3 = 1/Angstroms; R,D = Angstroms +# other quantities are unitless + +# This is the Ge reparameterization from the following paper: +# Optimized Tersoff empirical potential for germanene +# Mahdizadeh, Akhlamadi, J. Mol. Graph. Model. 72, 1-5 (2017) + +# format of a single entry (one or more lines): +# element 1, element 2, element 3, +# m, gamma, lambda3, c, d, costheta0, n, beta, lambda2, B, R, D, lambda1, A + +Ge Ge Ge 3.0 1.0 0.0 1.0643e5 15.2 -0.35 0.75627 5.017e-7 1.71 430.0 2.95 0.15 2.4451 1760.1 diff --git a/src/GRANULAR/fix_wall_gran_region.cpp b/src/GRANULAR/fix_wall_gran_region.cpp index d56715230..996708c1a 100644 --- a/src/GRANULAR/fix_wall_gran_region.cpp +++ b/src/GRANULAR/fix_wall_gran_region.cpp @@ -1,510 +1,510 @@ /* ---------------------------------------------------------------------- LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator http://lammps.sandia.gov, Sandia National Laboratories Steve Plimpton, sjplimp@sandia.gov Copyright (2003) Sandia Corporation. Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains certain rights in this software. This software is distributed under the GNU General Public License. See the README file in the top-level LAMMPS directory. ------------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- Contributing authors: Dan Bolintineanu (SNL) ------------------------------------------------------------------------- */ #include #include #include #include "fix_wall_gran_region.h" #include "region.h" #include "atom.h" #include "domain.h" #include "update.h" #include "force.h" #include "pair.h" #include "modify.h" #include "respa.h" #include "math_const.h" #include "memory.h" #include "error.h" using namespace LAMMPS_NS; using namespace FixConst; using namespace MathConst; // same as FixWallGran enum{HOOKE,HOOKE_HISTORY,HERTZ_HISTORY,BONDED_HISTORY}; #define BIG 1.0e20 /* ---------------------------------------------------------------------- */ FixWallGranRegion::FixWallGranRegion(LAMMPS *lmp, int narg, char **arg) : FixWallGran(lmp, narg, arg), region(NULL), region_style(NULL), ncontact(NULL), walls(NULL), shearmany(NULL), c2r(NULL) { restart_global = 1; motion_resetflag = 0; int iregion = domain->find_region(idregion); if (iregion == -1) error->all(FLERR,"Region ID for fix wall/gran/region does not exist"); region = domain->regions[iregion]; region_style = new char[strlen(region->style)+1]; strcpy(region_style,region->style); nregion = region->nregion; tmax = domain->regions[iregion]->tmax; c2r = new int[tmax]; // re-allocate atom-based arrays with nshear // do not register with Atom class, since parent class did that memory->destroy(shearone); shearone = NULL; ncontact = NULL; walls = NULL; shearmany = NULL; grow_arrays(atom->nmax); // initialize shear history as if particle is not touching region if (history) { int nlocal = atom->nlocal; for (int i = 0; i < nlocal; i++) ncontact[i] = 0; } } /* ---------------------------------------------------------------------- */ FixWallGranRegion::~FixWallGranRegion() { delete [] c2r; delete [] region_style; memory->destroy(ncontact); memory->destroy(walls); memory->destroy(shearmany); } /* ---------------------------------------------------------------------- */ void FixWallGranRegion::init() { FixWallGran::init(); int iregion = domain->find_region(idregion); if (iregion == -1) error->all(FLERR,"Region ID for fix wall/gran/region does not exist"); region = domain->regions[iregion]; // check if region properties changed between runs // reset if restart info was inconsistent if (strcmp(idregion,region->id) != 0 || strcmp(region_style,region->style) != 0 || nregion != region->nregion) { char str[256]; sprintf(str,"Region properties for region %s changed between runs, " "resetting its motion",idregion); error->warning(FLERR,str); region->reset_vel(); } if (motion_resetflag){ char str[256]; sprintf(str,"Region properties for region %s are inconsistent " "with restart file, resetting its motion",idregion); error->warning(FLERR,str); region->reset_vel(); } } /* ---------------------------------------------------------------------- */ void FixWallGranRegion::post_force(int vflag) { int i,m,nc,iwall; double rinv,fx,fy,fz,tooclose; double dx,dy,dz,rsq,meff; double xc[3],vwall[3]; // do not update shear history during setup shearupdate = 1; if (update->setupflag) shearupdate = 0; // if just reneighbored: // update rigid body masses for owned atoms if using FixRigid // body[i] = which body atom I is in, -1 if none // mass_body = mass of each rigid body if (neighbor->ago == 0 && fix_rigid) { int tmp; int *body = (int *) fix_rigid->extract("body",tmp); double *mass_body = (double *) fix_rigid->extract("masstotal",tmp); if (atom->nmax > nmax) { memory->destroy(mass_rigid); nmax = atom->nmax; memory->create(mass_rigid,nmax,"wall/gran:mass_rigid"); } int nlocal = atom->nlocal; for (i = 0; i < nlocal; i++) { if (body[i] >= 0) mass_rigid[i] = mass_body[body[i]]; else mass_rigid[i] = 0.0; } } int regiondynamic = region->dynamic_check(); if (!regiondynamic) vwall[0] = vwall[1] = vwall[2] = 0.0; double **x = atom->x; double **v = atom->v; double **f = atom->f; double **omega = atom->omega; double **torque = atom->torque; double *radius = atom->radius; double *rmass = atom->rmass; int *mask = atom->mask; int nlocal = atom->nlocal; // set current motion attributes of region // set_velocity() also updates prev to current step if (regiondynamic) { region->prematch(); region->set_velocity(); } for (i = 0; i < nlocal; i++) { if (mask[i] & groupbit) { if (!region->match(x[i][0],x[i][1],x[i][2])) continue; nc = region->surface(x[i][0],x[i][1],x[i][2],radius[i]); if (nc > tmax) error->one(FLERR,"Too many wall/gran/region contacts for one particle"); // shear history maintenance // update ncontact,walls,shear2many for particle I // to reflect new and persistent shear history values // also set c2r[] = indices into region->contact[] for each of N contacts // process zero or one contact here, otherwise invoke update_contacts() if (history) { if (nc == 0) { ncontact[i] = 0; continue; } if (nc == 1) { c2r[0] = 0; iwall = region->contact[0].iwall; if (ncontact[i] == 0) { ncontact[i] = 1; walls[i][0] = iwall; for (m = 0; m < sheardim; m++) shearmany[i][0][m] = 0.0; } else if (ncontact[i] > 1 || iwall != walls[i][0]) update_contacts(i,nc); } else update_contacts(i,nc); } // process current contacts for (int ic = 0; ic < nc; ic++) { // rsq = squared contact distance // xc = contact point rsq = region->contact[ic].r*region->contact[ic].r; dx = region->contact[ic].delx; dy = region->contact[ic].dely; dz = region->contact[ic].delz; if (regiondynamic) region->velocity_contact(vwall, x[i], ic); // meff = effective mass of sphere // if I is part of rigid body, use body mass meff = rmass[i]; if (fix_rigid && mass_rigid[i] > 0.0) meff = mass_rigid[i]; // invoke sphere/wall interaction if (pairstyle == HOOKE) hooke(rsq,dx,dy,dz,vwall,v[i],f[i], omega[i],torque[i],radius[i],meff); else if (pairstyle == HOOKE_HISTORY) hooke_history(rsq,dx,dy,dz,vwall,v[i],f[i], omega[i],torque[i],radius[i],meff, shearmany[i][c2r[ic]]); else if (pairstyle == HERTZ_HISTORY) hertz_history(rsq,dx,dy,dz,vwall,region->contact[ic].radius, v[i],f[i],omega[i],torque[i], radius[i],meff,shearmany[i][c2r[ic]]); else if (pairstyle == BONDED_HISTORY) bonded_history(rsq,dx,dy,dz,vwall,region->contact[ic].radius, v[i],f[i],omega[i],torque[i], radius[i],meff,shearmany[i][c2r[ic]]); } } } } /* ---------------------------------------------------------------------- update contact info in ncontact, walls, shear2many for particle I based on ncontacts[i] old contacts and N new contacts matched via their associated walls delete/zero shear history for broken/new contacts also set c2r[i] = index of Ith contact in region list of contacts ------------------------------------------------------------------------- */ void FixWallGranRegion::update_contacts(int i, int nc) { int j,m,iold,nold,ilast,inew,iadd,iwall; // loop over old contacts // if not in new contact list: // delete old contact by copying last contact over it iold = 0; while (iold < ncontact[i]) { for (m = 0; m < nc; m++) - if (region->contact[m].iwall = walls[i][iold]) break; + if (region->contact[m].iwall == walls[i][iold]) break; if (m < nc) { ilast = ncontact[i]-1; for (j = 0; j < sheardim; j++) shearmany[i][iold][j] = shearmany[i][ilast][j]; walls[i][iold] = walls[i][ilast]; ncontact[i]--; } else iold++; } // loop over new contacts // if not in newly compressed contact list of length nold: // add it with zeroed shear history // set all values in c2r nold = ncontact[i]; for (inew = 0; inew < nc; inew++) { iwall = region->contact[inew].iwall; for (m = 0; m < nold; m++) if (walls[i][m] == iwall) break; if (m < nold) c2r[m] = inew; else { iadd = ncontact[i]; c2r[iadd] = inew; for (j = 0; j < sheardim; j++) shearmany[i][iadd][j] = 0.0; walls[i][iadd] = iwall; ncontact[i]++; } } } /* ---------------------------------------------------------------------- memory usage of local atom-based arrays ------------------------------------------------------------------------- */ double FixWallGranRegion::memory_usage() { int nmax = atom->nmax; double bytes = 0.0; if (history) { // shear history bytes += nmax * sizeof(int); // ncontact bytes += nmax*tmax * sizeof(int); // walls bytes += nmax*tmax*sheardim * sizeof(double); // shearmany } if (fix_rigid) bytes += nmax * sizeof(int); // mass_rigid return bytes; } /* ---------------------------------------------------------------------- allocate local atom-based arrays ------------------------------------------------------------------------- */ void FixWallGranRegion::grow_arrays(int nmax) { if (history) { memory->grow(ncontact,nmax,"fix_wall_gran:ncontact"); memory->grow(walls,nmax,tmax,"fix_wall_gran:walls"); memory->grow(shearmany,nmax,tmax,sheardim,"fix_wall_gran:shearmany"); } } /* ---------------------------------------------------------------------- copy values within local atom-based arrays ------------------------------------------------------------------------- */ void FixWallGranRegion::copy_arrays(int i, int j, int delflag) { int m,n,iwall; if (!history) return; n = ncontact[i]; for (iwall = 0; iwall < n; iwall++) { walls[j][iwall] = walls[i][iwall]; for (m = 0; m < sheardim; m++) shearmany[j][iwall][m] = shearmany[i][iwall][m]; } ncontact[j] = ncontact[i]; } /* ---------------------------------------------------------------------- initialize one atom's array values, called when atom is created ------------------------------------------------------------------------- */ void FixWallGranRegion::set_arrays(int i) { if (!history) return; ncontact[i] = 0; } /* ---------------------------------------------------------------------- pack values in local atom-based arrays for exchange with another proc ------------------------------------------------------------------------- */ int FixWallGranRegion::pack_exchange(int i, double *buf) { int m; if (!history) return 0; int n = 0; int count = ncontact[i]; buf[n++] = ubuf(count).d; for (int iwall = 0; iwall < count; iwall++) { buf[n++] = ubuf(walls[i][iwall]).d; for (m = 0; m < sheardim; m++) buf[n++] = shearmany[i][iwall][m]; } return n; } /* ---------------------------------------------------------------------- unpack values into local atom-based arrays after exchange ------------------------------------------------------------------------- */ int FixWallGranRegion::unpack_exchange(int nlocal, double *buf) { int m; if (!history) return 0; int n = 0; int count = ncontact[nlocal] = (int) ubuf(buf[n++]).i; for (int iwall = 0; iwall < count; iwall++) { walls[nlocal][iwall] = (int) ubuf(buf[n++]).i; for (m = 0; m < sheardim; m++) shearmany[nlocal][iwall][m] = buf[n++]; } return n; } /* ---------------------------------------------------------------------- pack values in local atom-based arrays for restart file ------------------------------------------------------------------------- */ int FixWallGranRegion::pack_restart(int i, double *buf) { int m; if (!history) return 0; int n = 1; int count = ncontact[i]; buf[n++] = ubuf(count).d; for (int iwall = 0; iwall < count; iwall++) { buf[n++] = ubuf(walls[i][iwall]).d; for (m = 0; m < sheardim; m++) buf[n++] = shearmany[i][iwall][m]; } buf[0] = n; return n; } /* ---------------------------------------------------------------------- unpack values from atom->extra array to restart the fix ------------------------------------------------------------------------- */ void FixWallGranRegion::unpack_restart(int nlocal, int nth) { int k; if (!history) return; double **extra = atom->extra; // skip to Nth set of extra values int m = 0; for (int i = 0; i < nth; i++) m += static_cast (extra[nlocal][m]); m++; int count = ncontact[nlocal] = (int) ubuf(extra[nlocal][m++]).i; for (int iwall = 0; iwall < count; iwall++) { walls[nlocal][iwall] = (int) ubuf(extra[nlocal][m++]).i; for (k = 0; k < sheardim; k++) shearmany[nlocal][iwall][k] = extra[nlocal][m++]; } } /* ---------------------------------------------------------------------- maxsize of any atom's restart data ------------------------------------------------------------------------- */ int FixWallGranRegion::maxsize_restart() { if (!history) return 0; return 2 + tmax*(sheardim+1); } /* ---------------------------------------------------------------------- size of atom nlocal's restart data ------------------------------------------------------------------------- */ int FixWallGranRegion::size_restart(int nlocal) { if (!history) return 0; return 2 + ncontact[nlocal]*(sheardim+1); } /* ---------------------------------------------------------------------- pack entire state of Fix into one write ------------------------------------------------------------------------- */ void FixWallGranRegion::write_restart(FILE *fp) { if (comm->me) return; int len = 0; region->length_restart_string(len); fwrite(&len, sizeof(int),1,fp); region->write_restart(fp); } /* ---------------------------------------------------------------------- use state info from restart file to restart the Fix ------------------------------------------------------------------------- */ void FixWallGranRegion::restart(char *buf) { int n = 0; if (!region->restart(buf,n)) motion_resetflag = 1; } diff --git a/src/GRANULAR/pair_gran_hooke_history.cpp b/src/GRANULAR/pair_gran_hooke_history.cpp index edd12b4b2..ac5efb48d 100644 --- a/src/GRANULAR/pair_gran_hooke_history.cpp +++ b/src/GRANULAR/pair_gran_hooke_history.cpp @@ -1,801 +1,812 @@ /* ---------------------------------------------------------------------- LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator http://lammps.sandia.gov, Sandia National Laboratories Steve Plimpton, sjplimp@sandia.gov Copyright (2003) Sandia Corporation. Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains certain rights in this software. This software is distributed under the GNU General Public License. See the README file in the top-level LAMMPS directory. ------------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- Contributing authors: Leo Silbert (SNL), Gary Grest (SNL) ------------------------------------------------------------------------- */ #include #include #include #include #include "pair_gran_hooke_history.h" #include "atom.h" #include "atom_vec.h" #include "domain.h" #include "force.h" #include "update.h" #include "modify.h" #include "fix.h" #include "fix_shear_history.h" #include "comm.h" #include "neighbor.h" #include "neigh_list.h" #include "neigh_request.h" #include "memory.h" #include "error.h" using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ PairGranHookeHistory::PairGranHookeHistory(LAMMPS *lmp) : Pair(lmp) { single_enable = 1; no_virial_fdotr_compute = 1; history = 1; fix_history = NULL; single_extra = 10; svector = new double[10]; neighprev = 0; nmax = 0; mass_rigid = NULL; // set comm size needed by this Pair if used with fix rigid comm_forward = 1; } /* ---------------------------------------------------------------------- */ PairGranHookeHistory::~PairGranHookeHistory() { delete [] svector; if (fix_history) modify->delete_fix("SHEAR_HISTORY"); if (allocated) { memory->destroy(setflag); memory->destroy(cutsq); delete [] onerad_dynamic; delete [] onerad_frozen; delete [] maxrad_dynamic; delete [] maxrad_frozen; } memory->destroy(mass_rigid); } /* ---------------------------------------------------------------------- */ void PairGranHookeHistory::compute(int eflag, int vflag) { int i,j,ii,jj,inum,jnum; double xtmp,ytmp,ztmp,delx,dely,delz,fx,fy,fz; double radi,radj,radsum,rsq,r,rinv,rsqinv; double vr1,vr2,vr3,vnnr,vn1,vn2,vn3,vt1,vt2,vt3; double wr1,wr2,wr3; double vtr1,vtr2,vtr3,vrel; double mi,mj,meff,damp,ccel,tor1,tor2,tor3; double fn,fs,fs1,fs2,fs3; double shrmag,rsht; int *ilist,*jlist,*numneigh,**firstneigh; int *touch,**firsttouch; double *shear,*allshear,**firstshear; if (eflag || vflag) ev_setup(eflag,vflag); else evflag = vflag_fdotr = 0; int shearupdate = 1; if (update->setupflag) shearupdate = 0; // update rigid body info for owned & ghost atoms if using FixRigid masses // body[i] = which body atom I is in, -1 if none // mass_body = mass of each rigid body if (fix_rigid && neighbor->ago == 0) { int tmp; int *body = (int *) fix_rigid->extract("body",tmp); double *mass_body = (double *) fix_rigid->extract("masstotal",tmp); if (atom->nmax > nmax) { memory->destroy(mass_rigid); nmax = atom->nmax; memory->create(mass_rigid,nmax,"pair:mass_rigid"); } int nlocal = atom->nlocal; for (i = 0; i < nlocal; i++) if (body[i] >= 0) mass_rigid[i] = mass_body[body[i]]; else mass_rigid[i] = 0.0; comm->forward_comm_pair(this); } double **x = atom->x; double **v = atom->v; double **f = atom->f; double **omega = atom->omega; double **torque = atom->torque; double *radius = atom->radius; double *rmass = atom->rmass; int *mask = atom->mask; int nlocal = atom->nlocal; int newton_pair = force->newton_pair; inum = list->inum; ilist = list->ilist; numneigh = list->numneigh; firstneigh = list->firstneigh; firsttouch = listgranhistory->firstneigh; firstshear = listgranhistory->firstdouble; // loop over neighbors of my atoms for (ii = 0; ii < inum; ii++) { i = ilist[ii]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; radi = radius[i]; touch = firsttouch[i]; allshear = firstshear[i]; jlist = firstneigh[i]; jnum = numneigh[i]; for (jj = 0; jj < jnum; jj++) { j = jlist[jj]; j &= NEIGHMASK; delx = xtmp - x[j][0]; dely = ytmp - x[j][1]; delz = ztmp - x[j][2]; rsq = delx*delx + dely*dely + delz*delz; radj = radius[j]; radsum = radi + radj; if (rsq >= radsum*radsum) { // unset non-touching neighbors touch[jj] = 0; shear = &allshear[3*jj]; shear[0] = 0.0; shear[1] = 0.0; shear[2] = 0.0; } else { r = sqrt(rsq); rinv = 1.0/r; rsqinv = 1.0/rsq; // relative translational velocity vr1 = v[i][0] - v[j][0]; vr2 = v[i][1] - v[j][1]; vr3 = v[i][2] - v[j][2]; // normal component vnnr = vr1*delx + vr2*dely + vr3*delz; vn1 = delx*vnnr * rsqinv; vn2 = dely*vnnr * rsqinv; vn3 = delz*vnnr * rsqinv; // tangential component vt1 = vr1 - vn1; vt2 = vr2 - vn2; vt3 = vr3 - vn3; // relative rotational velocity wr1 = (radi*omega[i][0] + radj*omega[j][0]) * rinv; wr2 = (radi*omega[i][1] + radj*omega[j][1]) * rinv; wr3 = (radi*omega[i][2] + radj*omega[j][2]) * rinv; // meff = effective mass of pair of particles // if I or J part of rigid body, use body mass // if I or J is frozen, meff is other particle mi = rmass[i]; mj = rmass[j]; if (fix_rigid) { if (mass_rigid[i] > 0.0) mi = mass_rigid[i]; if (mass_rigid[j] > 0.0) mj = mass_rigid[j]; } meff = mi*mj / (mi+mj); if (mask[i] & freeze_group_bit) meff = mj; if (mask[j] & freeze_group_bit) meff = mi; // normal forces = Hookian contact + normal velocity damping damp = meff*gamman*vnnr*rsqinv; ccel = kn*(radsum-r)*rinv - damp; // relative velocities vtr1 = vt1 - (delz*wr2-dely*wr3); vtr2 = vt2 - (delx*wr3-delz*wr1); vtr3 = vt3 - (dely*wr1-delx*wr2); vrel = vtr1*vtr1 + vtr2*vtr2 + vtr3*vtr3; vrel = sqrt(vrel); // shear history effects touch[jj] = 1; shear = &allshear[3*jj]; if (shearupdate) { shear[0] += vtr1*dt; shear[1] += vtr2*dt; shear[2] += vtr3*dt; } shrmag = sqrt(shear[0]*shear[0] + shear[1]*shear[1] + shear[2]*shear[2]); // rotate shear displacements rsht = shear[0]*delx + shear[1]*dely + shear[2]*delz; rsht *= rsqinv; if (shearupdate) { shear[0] -= rsht*delx; shear[1] -= rsht*dely; shear[2] -= rsht*delz; } // tangential forces = shear + tangential velocity damping fs1 = - (kt*shear[0] + meff*gammat*vtr1); fs2 = - (kt*shear[1] + meff*gammat*vtr2); fs3 = - (kt*shear[2] + meff*gammat*vtr3); // rescale frictional displacements and forces if needed fs = sqrt(fs1*fs1 + fs2*fs2 + fs3*fs3); fn = xmu * fabs(ccel*r); if (fs > fn) { if (shrmag != 0.0) { shear[0] = (fn/fs) * (shear[0] + meff*gammat*vtr1/kt) - meff*gammat*vtr1/kt; shear[1] = (fn/fs) * (shear[1] + meff*gammat*vtr2/kt) - meff*gammat*vtr2/kt; shear[2] = (fn/fs) * (shear[2] + meff*gammat*vtr3/kt) - meff*gammat*vtr3/kt; fs1 *= fn/fs; fs2 *= fn/fs; fs3 *= fn/fs; } else fs1 = fs2 = fs3 = 0.0; } // forces & torques fx = delx*ccel + fs1; fy = dely*ccel + fs2; fz = delz*ccel + fs3; f[i][0] += fx; f[i][1] += fy; f[i][2] += fz; tor1 = rinv * (dely*fs3 - delz*fs2); tor2 = rinv * (delz*fs1 - delx*fs3); tor3 = rinv * (delx*fs2 - dely*fs1); torque[i][0] -= radi*tor1; torque[i][1] -= radi*tor2; torque[i][2] -= radi*tor3; if (newton_pair || j < nlocal) { f[j][0] -= fx; f[j][1] -= fy; f[j][2] -= fz; torque[j][0] -= radj*tor1; torque[j][1] -= radj*tor2; torque[j][2] -= radj*tor3; } if (evflag) ev_tally_xyz(i,j,nlocal,newton_pair, 0.0,0.0,fx,fy,fz,delx,dely,delz); } } } if (vflag_fdotr) virial_fdotr_compute(); } /* ---------------------------------------------------------------------- allocate all arrays ------------------------------------------------------------------------- */ void PairGranHookeHistory::allocate() { allocated = 1; int n = atom->ntypes; memory->create(setflag,n+1,n+1,"pair:setflag"); for (int i = 1; i <= n; i++) for (int j = i; j <= n; j++) setflag[i][j] = 0; memory->create(cutsq,n+1,n+1,"pair:cutsq"); onerad_dynamic = new double[n+1]; onerad_frozen = new double[n+1]; maxrad_dynamic = new double[n+1]; maxrad_frozen = new double[n+1]; } /* ---------------------------------------------------------------------- global settings ------------------------------------------------------------------------- */ void PairGranHookeHistory::settings(int narg, char **arg) { if (narg != 6) error->all(FLERR,"Illegal pair_style command"); kn = force->numeric(FLERR,arg[0]); if (strcmp(arg[1],"NULL") == 0) kt = kn * 2.0/7.0; else kt = force->numeric(FLERR,arg[1]); gamman = force->numeric(FLERR,arg[2]); if (strcmp(arg[3],"NULL") == 0) gammat = 0.5 * gamman; else gammat = force->numeric(FLERR,arg[3]); xmu = force->numeric(FLERR,arg[4]); dampflag = force->inumeric(FLERR,arg[5]); if (dampflag == 0) gammat = 0.0; if (kn < 0.0 || kt < 0.0 || gamman < 0.0 || gammat < 0.0 || xmu < 0.0 || xmu > 10000.0 || dampflag < 0 || dampflag > 1) error->all(FLERR,"Illegal pair_style command"); } /* ---------------------------------------------------------------------- set coeffs for one or more type pairs ------------------------------------------------------------------------- */ void PairGranHookeHistory::coeff(int narg, char **arg) { if (narg > 2) error->all(FLERR,"Incorrect args for pair coefficients"); if (!allocated) allocate(); int ilo,ihi,jlo,jhi; force->bounds(FLERR,arg[0],atom->ntypes,ilo,ihi); force->bounds(FLERR,arg[1],atom->ntypes,jlo,jhi); int count = 0; for (int i = ilo; i <= ihi; i++) { for (int j = MAX(jlo,i); j <= jhi; j++) { setflag[i][j] = 1; count++; } } if (count == 0) error->all(FLERR,"Incorrect args for pair coefficients"); } /* ---------------------------------------------------------------------- init specific to this pair style ------------------------------------------------------------------------- */ void PairGranHookeHistory::init_style() { int i; // error and warning checks if (!atom->radius_flag || !atom->rmass_flag) error->all(FLERR,"Pair granular requires atom attributes radius, rmass"); if (comm->ghost_velocity == 0) error->all(FLERR,"Pair granular requires ghost atoms store velocity"); // need a granular neigh list and optionally a granular history neigh list int irequest = neighbor->request(this,instance_me); neighbor->requests[irequest]->half = 0; neighbor->requests[irequest]->gran = 1; if (history) { irequest = neighbor->request(this,instance_me); neighbor->requests[irequest]->id = 1; neighbor->requests[irequest]->half = 0; neighbor->requests[irequest]->granhistory = 1; neighbor->requests[irequest]->dnum = 3; } dt = update->dt; // if shear history is stored: // if first init, create Fix needed for storing shear history if (history && fix_history == NULL) { char dnumstr[16]; sprintf(dnumstr,"%d",3); char **fixarg = new char*[4]; fixarg[0] = (char *) "SHEAR_HISTORY"; fixarg[1] = (char *) "all"; fixarg[2] = (char *) "SHEAR_HISTORY"; fixarg[3] = dnumstr; modify->add_fix(4,fixarg,1); delete [] fixarg; fix_history = (FixShearHistory *) modify->fix[modify->nfix-1]; fix_history->pair = this; neighbor->requests[irequest]->fix_history = fix_history; } // check for FixFreeze and set freeze_group_bit for (i = 0; i < modify->nfix; i++) if (strcmp(modify->fix[i]->style,"freeze") == 0) break; if (i < modify->nfix) freeze_group_bit = modify->fix[i]->groupbit; else freeze_group_bit = 0; // check for FixRigid so can extract rigid body masses fix_rigid = NULL; for (i = 0; i < modify->nfix; i++) if (modify->fix[i]->rigid_flag) break; if (i < modify->nfix) fix_rigid = modify->fix[i]; // check for FixPour and FixDeposit so can extract particle radii int ipour; for (ipour = 0; ipour < modify->nfix; ipour++) if (strcmp(modify->fix[ipour]->style,"pour") == 0) break; if (ipour == modify->nfix) ipour = -1; int idep; for (idep = 0; idep < modify->nfix; idep++) if (strcmp(modify->fix[idep]->style,"deposit") == 0) break; if (idep == modify->nfix) idep = -1; // set maxrad_dynamic and maxrad_frozen for each type // include future FixPour and FixDeposit particles as dynamic int itype; for (i = 1; i <= atom->ntypes; i++) { onerad_dynamic[i] = onerad_frozen[i] = 0.0; if (ipour >= 0) { itype = i; onerad_dynamic[i] = *((double *) modify->fix[ipour]->extract("radius",itype)); } if (idep >= 0) { itype = i; onerad_dynamic[i] = *((double *) modify->fix[idep]->extract("radius",itype)); } } double *radius = atom->radius; int *mask = atom->mask; int *type = atom->type; int nlocal = atom->nlocal; for (i = 0; i < nlocal; i++) if (mask[i] & freeze_group_bit) onerad_frozen[type[i]] = MAX(onerad_frozen[type[i]],radius[i]); else onerad_dynamic[type[i]] = MAX(onerad_dynamic[type[i]],radius[i]); MPI_Allreduce(&onerad_dynamic[1],&maxrad_dynamic[1],atom->ntypes, MPI_DOUBLE,MPI_MAX,world); MPI_Allreduce(&onerad_frozen[1],&maxrad_frozen[1],atom->ntypes, MPI_DOUBLE,MPI_MAX,world); // set fix which stores history info if (history) { int ifix = modify->find_fix("SHEAR_HISTORY"); if (ifix < 0) error->all(FLERR,"Could not find pair fix ID"); fix_history = (FixShearHistory *) modify->fix[ifix]; } } /* ---------------------------------------------------------------------- neighbor callback to inform pair style of neighbor list to use optional granular history list ------------------------------------------------------------------------- */ void PairGranHookeHistory::init_list(int id, NeighList *ptr) { if (id == 0) list = ptr; else if (id == 1) listgranhistory = ptr; } /* ---------------------------------------------------------------------- init for one type pair i,j and corresponding j,i ------------------------------------------------------------------------- */ double PairGranHookeHistory::init_one(int i, int j) { if (!allocated) allocate(); // cutoff = sum of max I,J radii for // dynamic/dynamic & dynamic/frozen interactions, but not frozen/frozen double cutoff = maxrad_dynamic[i]+maxrad_dynamic[j]; cutoff = MAX(cutoff,maxrad_frozen[i]+maxrad_dynamic[j]); cutoff = MAX(cutoff,maxrad_dynamic[i]+maxrad_frozen[j]); return cutoff; } /* ---------------------------------------------------------------------- proc 0 writes to restart file ------------------------------------------------------------------------- */ void PairGranHookeHistory::write_restart(FILE *fp) { write_restart_settings(fp); int i,j; for (i = 1; i <= atom->ntypes; i++) for (j = i; j <= atom->ntypes; j++) fwrite(&setflag[i][j],sizeof(int),1,fp); } /* ---------------------------------------------------------------------- proc 0 reads from restart file, bcasts ------------------------------------------------------------------------- */ void PairGranHookeHistory::read_restart(FILE *fp) { read_restart_settings(fp); allocate(); int i,j; int me = comm->me; for (i = 1; i <= atom->ntypes; i++) for (j = i; j <= atom->ntypes; j++) { if (me == 0) fread(&setflag[i][j],sizeof(int),1,fp); MPI_Bcast(&setflag[i][j],1,MPI_INT,0,world); } } /* ---------------------------------------------------------------------- proc 0 writes to restart file ------------------------------------------------------------------------- */ void PairGranHookeHistory::write_restart_settings(FILE *fp) { fwrite(&kn,sizeof(double),1,fp); fwrite(&kt,sizeof(double),1,fp); fwrite(&gamman,sizeof(double),1,fp); fwrite(&gammat,sizeof(double),1,fp); fwrite(&xmu,sizeof(double),1,fp); fwrite(&dampflag,sizeof(int),1,fp); } /* ---------------------------------------------------------------------- proc 0 reads from restart file, bcasts ------------------------------------------------------------------------- */ void PairGranHookeHistory::read_restart_settings(FILE *fp) { if (comm->me == 0) { fread(&kn,sizeof(double),1,fp); fread(&kt,sizeof(double),1,fp); fread(&gamman,sizeof(double),1,fp); fread(&gammat,sizeof(double),1,fp); fread(&xmu,sizeof(double),1,fp); fread(&dampflag,sizeof(int),1,fp); } MPI_Bcast(&kn,1,MPI_DOUBLE,0,world); MPI_Bcast(&kt,1,MPI_DOUBLE,0,world); MPI_Bcast(&gamman,1,MPI_DOUBLE,0,world); MPI_Bcast(&gammat,1,MPI_DOUBLE,0,world); MPI_Bcast(&xmu,1,MPI_DOUBLE,0,world); MPI_Bcast(&dampflag,1,MPI_INT,0,world); } /* ---------------------------------------------------------------------- */ void PairGranHookeHistory::reset_dt() { dt = update->dt; } /* ---------------------------------------------------------------------- */ double PairGranHookeHistory::single(int i, int j, int itype, int jtype, double rsq, double factor_coul, double factor_lj, double &fforce) { double radi,radj,radsum; double r,rinv,rsqinv,delx,dely,delz; double vr1,vr2,vr3,vnnr,vn1,vn2,vn3,vt1,vt2,vt3,wr1,wr2,wr3; double mi,mj,meff,damp,ccel; double vtr1,vtr2,vtr3,vrel,shrmag,rsht; double fs1,fs2,fs3,fs,fn; double *radius = atom->radius; radi = radius[i]; radj = radius[j]; radsum = radi + radj; if (rsq >= radsum*radsum) { fforce = 0.0; for (int m = 0; m < single_extra; m++) svector[m] = 0.0; return 0.0; } r = sqrt(rsq); rinv = 1.0/r; rsqinv = 1.0/rsq; // relative translational velocity double **v = atom->v; vr1 = v[i][0] - v[j][0]; vr2 = v[i][1] - v[j][1]; vr3 = v[i][2] - v[j][2]; // normal component double **x = atom->x; delx = x[i][0] - x[j][0]; dely = x[i][1] - x[j][1]; delz = x[i][2] - x[j][2]; vnnr = vr1*delx + vr2*dely + vr3*delz; vn1 = delx*vnnr * rsqinv; vn2 = dely*vnnr * rsqinv; vn3 = delz*vnnr * rsqinv; // tangential component vt1 = vr1 - vn1; vt2 = vr2 - vn2; vt3 = vr3 - vn3; // relative rotational velocity double **omega = atom->omega; wr1 = (radi*omega[i][0] + radj*omega[j][0]) * rinv; wr2 = (radi*omega[i][1] + radj*omega[j][1]) * rinv; wr3 = (radi*omega[i][2] + radj*omega[j][2]) * rinv; // meff = effective mass of pair of particles // if I or J part of rigid body, use body mass // if I or J is frozen, meff is other particle double *rmass = atom->rmass; int *mask = atom->mask; mi = rmass[i]; mj = rmass[j]; if (fix_rigid) { // NOTE: insure mass_rigid is current for owned+ghost atoms? if (mass_rigid[i] > 0.0) mi = mass_rigid[i]; if (mass_rigid[j] > 0.0) mj = mass_rigid[j]; } meff = mi*mj / (mi+mj); if (mask[i] & freeze_group_bit) meff = mj; if (mask[j] & freeze_group_bit) meff = mi; // normal forces = Hookian contact + normal velocity damping damp = meff*gamman*vnnr*rsqinv; ccel = kn*(radsum-r)*rinv - damp; // relative velocities vtr1 = vt1 - (delz*wr2-dely*wr3); vtr2 = vt2 - (delx*wr3-delz*wr1); vtr3 = vt3 - (dely*wr1-delx*wr2); vrel = vtr1*vtr1 + vtr2*vtr2 + vtr3*vtr3; vrel = sqrt(vrel); // shear history effects // neighprev = index of found neigh on previous call // search entire jnum list of neighbors of I for neighbor J // start from neighprev, since will typically be next neighbor // reset neighprev to 0 as necessary int jnum = list->numneigh[i]; int *jlist = list->firstneigh[i]; double *allshear = list->listgranhistory->firstdouble[i]; for (int jj = 0; jj < jnum; jj++) { neighprev++; if (neighprev >= jnum) neighprev = 0; if (jlist[neighprev] == j) break; } double *shear = &allshear[3*neighprev]; shrmag = sqrt(shear[0]*shear[0] + shear[1]*shear[1] + shear[2]*shear[2]); // rotate shear displacements rsht = shear[0]*delx + shear[1]*dely + shear[2]*delz; rsht *= rsqinv; // tangential forces = shear + tangential velocity damping fs1 = - (kt*shear[0] + meff*gammat*vtr1); fs2 = - (kt*shear[1] + meff*gammat*vtr2); fs3 = - (kt*shear[2] + meff*gammat*vtr3); // rescale frictional displacements and forces if needed fs = sqrt(fs1*fs1 + fs2*fs2 + fs3*fs3); fn = xmu * fabs(ccel*r); if (fs > fn) { if (shrmag != 0.0) { fs1 *= fn/fs; fs2 *= fn/fs; fs3 *= fn/fs; fs *= fn/fs; } else fs1 = fs2 = fs3 = fs = 0.0; } // set force and return no energy fforce = ccel; // set single_extra quantities svector[0] = fs1; svector[1] = fs2; svector[2] = fs3; svector[3] = fs; svector[4] = vn1; svector[5] = vn2; svector[6] = vn3; svector[7] = vt1; svector[8] = vt2; svector[9] = vt3; return 0.0; } /* ---------------------------------------------------------------------- */ int PairGranHookeHistory::pack_forward_comm(int n, int *list, double *buf, int pbc_flag, int *pbc) { int i,j,m; m = 0; for (i = 0; i < n; i++) { j = list[i]; buf[m++] = mass_rigid[j]; } return m; } /* ---------------------------------------------------------------------- */ void PairGranHookeHistory::unpack_forward_comm(int n, int first, double *buf) { int i,m,last; m = 0; last = first + n; for (i = first; i < last; i++) mass_rigid[i] = buf[m++]; } /* ---------------------------------------------------------------------- memory usage of local atom-based arrays ------------------------------------------------------------------------- */ double PairGranHookeHistory::memory_usage() { double bytes = nmax * sizeof(double); return bytes; } + +/* ---------------------------------------------------------------------- + return ptr to FixShearHistory class + called by Neighbor when setting up neighbor lists +------------------------------------------------------------------------- */ + +void *PairGranHookeHistory::extract(const char *str, int &dim) +{ + if (strcmp(str,"history") == 0) return (void *) fix_history; + return NULL; +} diff --git a/src/GRANULAR/pair_gran_hooke_history.h b/src/GRANULAR/pair_gran_hooke_history.h index ecd52f353..afeab9341 100644 --- a/src/GRANULAR/pair_gran_hooke_history.h +++ b/src/GRANULAR/pair_gran_hooke_history.h @@ -1,104 +1,105 @@ /* -*- 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. ------------------------------------------------------------------------- */ #ifdef PAIR_CLASS PairStyle(gran/hooke/history,PairGranHookeHistory) #else #ifndef LMP_PAIR_GRAN_HOOKE_HISTORY_H #define LMP_PAIR_GRAN_HOOKE_HISTORY_H #include "pair.h" namespace LAMMPS_NS { class PairGranHookeHistory : public Pair { public: PairGranHookeHistory(class LAMMPS *); virtual ~PairGranHookeHistory(); virtual void compute(int, int); virtual void settings(int, char **); void coeff(int, char **); void init_style(); void init_list(int, class NeighList *); double init_one(int, int); void write_restart(FILE *); void read_restart(FILE *); void write_restart_settings(FILE *); void read_restart_settings(FILE *); void reset_dt(); virtual double single(int, int, int, int, double, double, double, double &); int pack_forward_comm(int, int *, double *, int, int *); void unpack_forward_comm(int, int, double *); double memory_usage(); + void *extract(const char *, int &); protected: double kn,kt,gamman,gammat,xmu; int dampflag; double dt; int freeze_group_bit; int history; int neighprev; double *onerad_dynamic,*onerad_frozen; double *maxrad_dynamic,*maxrad_frozen; class FixShearHistory *fix_history; // storage of rigid body masses for use in granular interactions class Fix *fix_rigid; // ptr to rigid body fix, NULL if none double *mass_rigid; // rigid mass for owned+ghost atoms int nmax; // allocated size of mass_rigid void allocate(); }; } #endif #endif /* ERROR/WARNING messages: E: Illegal ... command Self-explanatory. Check the input script syntax and compare to the documentation for the command. You can use -echo screen as a command-line option when running LAMMPS to see the offending line. E: Incorrect args for pair coefficients Self-explanatory. Check the input script or data file. E: Pair granular requires atom attributes radius, rmass The atom style defined does not have these attributes. E: Pair granular requires ghost atoms store velocity Use the comm_modify vel yes command to enable this. E: Pair granular with shear history requires newton pair off This is a current restriction of the implementation of pair granular styles with history. E: Could not find pair fix ID A fix is created internally by the pair style to store shear history information. You cannot delete it. */ diff --git a/src/USER-MISC/dihedral_nharmonic.cpp b/src/USER-MISC/dihedral_nharmonic.cpp index 87f7b0775..8958340bb 100644 --- a/src/USER-MISC/dihedral_nharmonic.cpp +++ b/src/USER-MISC/dihedral_nharmonic.cpp @@ -1,353 +1,355 @@ /* ---------------------------------------------------------------------- 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: Loukas D. Peristeras (Scienomics SARL) [ based on dihedral_multi_harmonic.cpp Mathias Puetz (SNL) and friends ] ------------------------------------------------------------------------- */ #include #include #include "dihedral_nharmonic.h" #include "atom.h" #include "neighbor.h" #include "domain.h" #include "comm.h" #include "force.h" #include "update.h" #include "memory.h" #include "error.h" using namespace LAMMPS_NS; #define TOLERANCE 0.05 #define SMALL 0.001 /* ---------------------------------------------------------------------- */ DihedralNHarmonic::DihedralNHarmonic(LAMMPS *lmp) : Dihedral(lmp) { writedata = 1; } /* ---------------------------------------------------------------------- */ DihedralNHarmonic::~DihedralNHarmonic() { if (allocated) { memory->destroy(setflag); for (int i = 1; i <= atom->ndihedraltypes; i++) delete [] a[i]; delete [] a; delete [] nterms; } } /* ---------------------------------------------------------------------- */ void DihedralNHarmonic::compute(int eflag, int vflag) { int i1,i2,i3,i4,n,type; double vb1x,vb1y,vb1z,vb2x,vb2y,vb2z,vb3x,vb3y,vb3z,vb2xm,vb2ym,vb2zm; double edihedral,f1[3],f2[3],f3[3],f4[3]; double sb1,sb2,sb3,rb1,rb3,c0,b1mag2,b1mag,b2mag2; double b2mag,b3mag2,b3mag,ctmp,c_,r12c1,c1mag,r12c2; double c2mag,sc1,sc2,s1,s12,c,p,pd,a11,a22; double a33,a12,a13,a23,sx2,sy2,sz2; double s2,sin2; edihedral = 0.0; if (eflag || vflag) ev_setup(eflag,vflag); else evflag = 0; double **x = atom->x; double **f = atom->f; int **dihedrallist = neighbor->dihedrallist; int ndihedrallist = neighbor->ndihedrallist; int nlocal = atom->nlocal; int newton_bond = force->newton_bond; for (n = 0; n < ndihedrallist; n++) { i1 = dihedrallist[n][0]; i2 = dihedrallist[n][1]; i3 = dihedrallist[n][2]; i4 = dihedrallist[n][3]; type = dihedrallist[n][4]; // 1st bond vb1x = x[i1][0] - x[i2][0]; vb1y = x[i1][1] - x[i2][1]; vb1z = x[i1][2] - x[i2][2]; // 2nd bond vb2x = x[i3][0] - x[i2][0]; vb2y = x[i3][1] - x[i2][1]; vb2z = x[i3][2] - x[i2][2]; vb2xm = -vb2x; vb2ym = -vb2y; vb2zm = -vb2z; // 3rd bond vb3x = x[i4][0] - x[i3][0]; vb3y = x[i4][1] - x[i3][1]; vb3z = x[i4][2] - x[i3][2]; // c0 calculation sb1 = 1.0 / (vb1x*vb1x + vb1y*vb1y + vb1z*vb1z); sb2 = 1.0 / (vb2x*vb2x + vb2y*vb2y + vb2z*vb2z); sb3 = 1.0 / (vb3x*vb3x + vb3y*vb3y + vb3z*vb3z); rb1 = sqrt(sb1); rb3 = sqrt(sb3); c0 = (vb1x*vb3x + vb1y*vb3y + vb1z*vb3z) * rb1*rb3; // 1st and 2nd angle b1mag2 = vb1x*vb1x + vb1y*vb1y + vb1z*vb1z; b1mag = sqrt(b1mag2); b2mag2 = vb2x*vb2x + vb2y*vb2y + vb2z*vb2z; b2mag = sqrt(b2mag2); b3mag2 = vb3x*vb3x + vb3y*vb3y + vb3z*vb3z; b3mag = sqrt(b3mag2); ctmp = vb1x*vb2x + vb1y*vb2y + vb1z*vb2z; r12c1 = 1.0 / (b1mag*b2mag); c1mag = ctmp * r12c1; ctmp = vb2xm*vb3x + vb2ym*vb3y + vb2zm*vb3z; r12c2 = 1.0 / (b2mag*b3mag); c2mag = ctmp * r12c2; // cos and sin of 2 angles and final c sin2 = MAX(1.0 - c1mag*c1mag,0.0); sc1 = sqrt(sin2); if (sc1 < SMALL) sc1 = SMALL; sc1 = 1.0/sc1; sin2 = MAX(1.0 - c2mag*c2mag,0.0); sc2 = sqrt(sin2); if (sc2 < SMALL) sc2 = SMALL; sc2 = 1.0/sc2; s1 = sc1 * sc1; s2 = sc2 * sc2; s12 = sc1 * sc2; c = (c0 + c1mag*c2mag) * s12; // error check if (c > 1.0 + TOLERANCE || c < (-1.0 - TOLERANCE)) { int me; MPI_Comm_rank(world,&me); if (screen) { char str[128]; sprintf(str,"Dihedral problem: %d " BIGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT " " TAGINT_FORMAT, me,update->ntimestep, atom->tag[i1],atom->tag[i2],atom->tag[i3],atom->tag[i4]); error->warning(FLERR,str,0); fprintf(screen," 1st atom: %d %g %g %g\n", me,x[i1][0],x[i1][1],x[i1][2]); fprintf(screen," 2nd atom: %d %g %g %g\n", me,x[i2][0],x[i2][1],x[i2][2]); fprintf(screen," 3rd atom: %d %g %g %g\n", me,x[i3][0],x[i3][1],x[i3][2]); fprintf(screen," 4th atom: %d %g %g %g\n", me,x[i4][0],x[i4][1],x[i4][2]); } } if (c > 1.0) c = 1.0; if (c < -1.0) c = -1.0; // force & energy // p = sum (i=1,n) a_i * c**(i-1) // pd = dp/dc + c_ = c; p = this->a[type][0]; pd = this->a[type][1]; for (int i = 1; i < nterms[type]-1; i++) { p += c_ * this->a[type][i]; pd += c_ * static_cast(i+1) * this->a[type][i+1]; c_ *= c; } p += c_ * this->a[type][nterms[type]-1]; if (eflag) edihedral = p; c = c * pd; s12 = s12 * pd; a11 = c*sb1*s1; a22 = -sb2 * (2.0*c0*s12 - c*(s1+s2)); a33 = c*sb3*s2; a12 = -r12c1*(c1mag*c*s1 + c2mag*s12); a13 = -rb1*rb3*s12; a23 = r12c2*(c2mag*c*s2 + c1mag*s12); sx2 = a12*vb1x + a22*vb2x + a23*vb3x; sy2 = a12*vb1y + a22*vb2y + a23*vb3y; sz2 = a12*vb1z + a22*vb2z + a23*vb3z; f1[0] = a11*vb1x + a12*vb2x + a13*vb3x; f1[1] = a11*vb1y + a12*vb2y + a13*vb3y; f1[2] = a11*vb1z + a12*vb2z + a13*vb3z; f2[0] = -sx2 - f1[0]; f2[1] = -sy2 - f1[1]; f2[2] = -sz2 - f1[2]; f4[0] = a13*vb1x + a23*vb2x + a33*vb3x; f4[1] = a13*vb1y + a23*vb2y + a33*vb3y; f4[2] = a13*vb1z + a23*vb2z + a33*vb3z; f3[0] = sx2 - f4[0]; f3[1] = sy2 - f4[1]; f3[2] = sz2 - f4[2]; // apply force to each of 4 atoms if (newton_bond || i1 < nlocal) { f[i1][0] += f1[0]; f[i1][1] += f1[1]; f[i1][2] += f1[2]; } if (newton_bond || i2 < nlocal) { f[i2][0] += f2[0]; f[i2][1] += f2[1]; f[i2][2] += f2[2]; } if (newton_bond || i3 < nlocal) { f[i3][0] += f3[0]; f[i3][1] += f3[1]; f[i3][2] += f3[2]; } if (newton_bond || i4 < nlocal) { f[i4][0] += f4[0]; f[i4][1] += f4[1]; f[i4][2] += f4[2]; } if (evflag) ev_tally(i1,i2,i3,i4,nlocal,newton_bond,edihedral,f1,f3,f4, vb1x,vb1y,vb1z,vb2x,vb2y,vb2z,vb3x,vb3y,vb3z); } } /* ---------------------------------------------------------------------- */ void DihedralNHarmonic::allocate() { allocated = 1; int n = atom->ndihedraltypes; memory->create(nterms,n+1,"dihedral:nt"); a = new double * [n+1]; memory->create(setflag,n+1,"dihedral:setflag"); for (int i = 1; i <= n; i++) setflag[i] = 0; } /* ---------------------------------------------------------------------- set coeffs for one type ------------------------------------------------------------------------- */ void DihedralNHarmonic::coeff(int narg, char **arg) { if (narg < 4 ) error->all(FLERR,"Incorrect args for dihedral coefficients"); int n = force->inumeric(FLERR,arg[1]); - if (narg != n + 2 ) error->all(FLERR,"Incorrect args for dihedral coefficients"); + if (narg != n + 2) + error->all(FLERR,"Incorrect args for dihedral coefficients"); if (!allocated) allocate(); int ilo,ihi; force->bounds(FLERR,arg[0],atom->ndihedraltypes,ilo,ihi); int count = 0; for (int i = ilo; i <= ihi; i++) { a[i] = new double [n]; nterms[i] = n; for (int j = 0; j < n; j++ ) { a[i][j] = force->numeric(FLERR,arg[2+j]); setflag[i] = 1; } count++; } if (count == 0) error->all(FLERR,"Incorrect args for dihedral coefficients"); } /* ---------------------------------------------------------------------- proc 0 writes out coeffs to restart file ------------------------------------------------------------------------- */ void DihedralNHarmonic::write_restart(FILE *fp) { fwrite(&nterms[1],sizeof(int),atom->ndihedraltypes,fp); for(int i = 1; i <= atom->ndihedraltypes; i++) fwrite(a[i],sizeof(double),nterms[i],fp); } /* ---------------------------------------------------------------------- proc 0 reads coeffs from restart file, bcasts them ------------------------------------------------------------------------- */ void DihedralNHarmonic::read_restart(FILE *fp) { allocate(); if (comm->me == 0) fread(&nterms[1],sizeof(int),atom->ndihedraltypes,fp); MPI_Bcast(&nterms[1],atom->ndihedraltypes,MPI_INT,0,world); // allocate for(int i = 1; i <= atom->ndihedraltypes; i++) a[i] = new double [nterms[i]]; if (comm->me == 0) { for(int i = 1; i <= atom->ndihedraltypes; i++) fread(a[i],sizeof(double),nterms[i],fp); } for (int i = 1; i <= atom->ndihedraltypes; i++ ) MPI_Bcast(a[i],nterms[i],MPI_DOUBLE,0,world); for (int i = 1; i <= atom->ndihedraltypes; i++) setflag[i] = 1; } /* ---------------------------------------------------------------------- proc 0 writes to data file ------------------------------------------------------------------------- */ void DihedralNHarmonic::write_data(FILE *fp) { for (int i = 1; i <= atom->ndihedraltypes; i++) { fprintf(fp, "%d %d", i, nterms[i]); for (int j = 0; j < nterms[i]; j++ ) fprintf(fp, " %g", a[i][j]); fprintf(fp, "\n"); } } diff --git a/src/compute_group_group.cpp b/src/compute_group_group.cpp index 4b5e6ad05..3d76c34da 100644 --- a/src/compute_group_group.cpp +++ b/src/compute_group_group.cpp @@ -1,402 +1,430 @@ /* ---------------------------------------------------------------------- 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: Naveen Michaud-Agrawal (Johns Hopkins U) K-space terms added by Stan Moore (BYU) ------------------------------------------------------------------------- */ #include #include #include "compute_group_group.h" #include "atom.h" #include "update.h" #include "force.h" #include "pair.h" #include "neighbor.h" #include "neigh_request.h" #include "neigh_list.h" #include "group.h" #include "kspace.h" #include "error.h" #include #include "comm.h" #include "domain.h" #include "math_const.h" using namespace LAMMPS_NS; using namespace MathConst; #define SMALL 0.00001 +enum{OFF,INTER,INTRA}; + /* ---------------------------------------------------------------------- */ ComputeGroupGroup::ComputeGroupGroup(LAMMPS *lmp, int narg, char **arg) : Compute(lmp, narg, arg), group2(NULL) { if (narg < 4) error->all(FLERR,"Illegal compute group/group command"); scalar_flag = vector_flag = 1; size_vector = 3; extscalar = 1; extvector = 1; int n = strlen(arg[3]) + 1; group2 = new char[n]; strcpy(group2,arg[3]); jgroup = group->find(group2); if (jgroup == -1) error->all(FLERR,"Compute group/group group ID does not exist"); jgroupbit = group->bitmask[jgroup]; pairflag = 1; kspaceflag = 0; boundaryflag = 1; + molflag = OFF; int iarg = 4; while (iarg < narg) { if (strcmp(arg[iarg],"pair") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal compute group/group command"); if (strcmp(arg[iarg+1],"yes") == 0) pairflag = 1; else if (strcmp(arg[iarg+1],"no") == 0) pairflag = 0; else error->all(FLERR,"Illegal compute group/group command"); iarg += 2; } else if (strcmp(arg[iarg],"kspace") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal compute group/group command"); if (strcmp(arg[iarg+1],"yes") == 0) kspaceflag = 1; else if (strcmp(arg[iarg+1],"no") == 0) kspaceflag = 0; else error->all(FLERR,"Illegal compute group/group command"); iarg += 2; } else if (strcmp(arg[iarg],"boundary") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal compute group/group command"); if (strcmp(arg[iarg+1],"yes") == 0) boundaryflag = 1; else if (strcmp(arg[iarg+1],"no") == 0) boundaryflag = 0; else error->all(FLERR,"Illegal compute group/group command"); iarg += 2; + } else if (strcmp(arg[iarg],"molecule") == 0) { + if (iarg+2 > narg) + error->all(FLERR,"Illegal compute group/group command"); + if (strcmp(arg[iarg+1],"off") == 0) molflag = OFF; + else if (strcmp(arg[iarg+1],"inter") == 0) molflag = INTER; + else if (strcmp(arg[iarg+1],"intra") == 0) molflag = INTRA; + else error->all(FLERR,"Illegal compute group/group command"); + if (molflag != OFF && atom->molecule_flag == 0) + error->all(FLERR,"Compute group/group molecule requires molecule IDs"); + iarg += 2; } else error->all(FLERR,"Illegal compute group/group command"); } vector = new double[3]; } /* ---------------------------------------------------------------------- */ ComputeGroupGroup::~ComputeGroupGroup() { delete [] group2; delete [] vector; } /* ---------------------------------------------------------------------- */ void ComputeGroupGroup::init() { // if non-hybrid, then error if single_enable = 0 // if hybrid, let hybrid determine if sub-style sets single_enable = 0 if (pairflag && force->pair == NULL) error->all(FLERR,"No pair style defined for compute group/group"); if (force->pair_match("hybrid",0) == NULL && force->pair->single_enable == 0) error->all(FLERR,"Pair style does not support compute group/group"); // error if Kspace style does not compute group/group interactions if (kspaceflag && force->kspace == NULL) error->all(FLERR,"No Kspace style defined for compute group/group"); if (kspaceflag && force->kspace->group_group_enable == 0) error->all(FLERR,"Kspace style does not support compute group/group"); if (pairflag) { pair = force->pair; cutsq = force->pair->cutsq; } else pair = NULL; if (kspaceflag) kspace = force->kspace; else kspace = NULL; // compute Kspace correction terms if (kspaceflag) { kspace_correction(); if (fabs(e_correction) > SMALL && comm->me == 0) { char str[128]; sprintf(str,"Both groups in compute group/group have a net charge; " "the Kspace boundary correction to energy will be non-zero"); error->warning(FLERR,str); } } // recheck that group 2 has not been deleted jgroup = group->find(group2); if (jgroup == -1) error->all(FLERR,"Compute group/group group ID does not exist"); jgroupbit = group->bitmask[jgroup]; // need an occasional half neighbor list if (pairflag) { int irequest = neighbor->request(this,instance_me); neighbor->requests[irequest]->pair = 0; neighbor->requests[irequest]->compute = 1; neighbor->requests[irequest]->occasional = 1; } } /* ---------------------------------------------------------------------- */ void ComputeGroupGroup::init_list(int id, NeighList *ptr) { list = ptr; } /* ---------------------------------------------------------------------- */ double ComputeGroupGroup::compute_scalar() { invoked_scalar = invoked_vector = update->ntimestep; scalar = 0.0; vector[0] = vector[1] = vector[2] = 0.0; if (pairflag) pair_contribution(); if (kspaceflag) kspace_contribution(); return scalar; } /* ---------------------------------------------------------------------- */ void ComputeGroupGroup::compute_vector() { invoked_scalar = invoked_vector = update->ntimestep; scalar = 0.0; vector[0] = vector[1] = vector[2] = 0.0; if (pairflag) pair_contribution(); if (kspaceflag) kspace_contribution(); } /* ---------------------------------------------------------------------- */ void ComputeGroupGroup::pair_contribution() { int i,j,ii,jj,inum,jnum,itype,jtype; double xtmp,ytmp,ztmp,delx,dely,delz; double rsq,eng,fpair,factor_coul,factor_lj; int *ilist,*jlist,*numneigh,**firstneigh; double **x = atom->x; + int *molecule = atom->molecule; int *type = atom->type; int *mask = atom->mask; int nlocal = atom->nlocal; double *special_coul = force->special_coul; double *special_lj = force->special_lj; int newton_pair = force->newton_pair; // invoke half neighbor list (will copy or build if necessary) neighbor->build_one(list); inum = list->inum; ilist = list->ilist; numneigh = list->numneigh; firstneigh = list->firstneigh; // loop over neighbors of my atoms // skip if I,J are not in 2 groups double one[4]; one[0] = one[1] = one[2] = one[3] = 0.0; for (ii = 0; ii < inum; ii++) { i = ilist[ii]; // skip if atom I is not in either group if (!(mask[i] & groupbit || mask[i] & jgroupbit)) continue; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; itype = type[i]; jlist = firstneigh[i]; jnum = numneigh[i]; for (jj = 0; jj < jnum; jj++) { j = jlist[jj]; factor_lj = special_lj[sbmask(j)]; factor_coul = special_coul[sbmask(j)]; j &= NEIGHMASK; - if (!(mask[j] & groupbit || mask[j] & jgroupbit)) continue; // skip if atom J is not in either group + // skip if atom J is not in either group + + if (!(mask[j] & groupbit || mask[j] & jgroupbit)) continue; + + // skip if atoms I,J are only in the same group int ij_flag = 0; int ji_flag = 0; if (mask[i] & groupbit && mask[j] & jgroupbit) ij_flag = 1; if (mask[j] & groupbit && mask[i] & jgroupbit) ji_flag = 1; - if (!ij_flag && !ji_flag) continue; // skip if atoms I,J are only in the same group + if (!ij_flag && !ji_flag) continue; + + // skip if molecule IDs of atoms I,J do not satisfy molflag setting + + if (molflag != OFF) { + if (molflag == INTER) { + if (molecule[i] == molecule[j]) continue; + } else { + if (molecule[i] != molecule[j]) continue; + } + } delx = xtmp - x[j][0]; dely = ytmp - x[j][1]; delz = ztmp - x[j][2]; rsq = delx*delx + dely*dely + delz*delz; jtype = type[j]; if (rsq < cutsq[itype][jtype]) { eng = pair->single(i,j,itype,jtype,rsq,factor_coul,factor_lj,fpair); // energy only computed once so tally full amount // force tally is jgroup acting on igroup if (newton_pair || j < nlocal) { one[0] += eng; if (ij_flag) { one[1] += delx*fpair; one[2] += dely*fpair; one[3] += delz*fpair; } if (ji_flag) { one[1] -= delx*fpair; one[2] -= dely*fpair; one[3] -= delz*fpair; } // energy computed twice so tally half amount // only tally force if I own igroup atom } else { one[0] += 0.5*eng; if (ij_flag) { one[1] += delx*fpair; one[2] += dely*fpair; one[3] += delz*fpair; } } } } } double all[4]; MPI_Allreduce(one,all,4,MPI_DOUBLE,MPI_SUM,world); scalar += all[0]; vector[0] += all[1]; vector[1] += all[2]; vector[2] += all[3]; } /* ---------------------------------------------------------------------- */ void ComputeGroupGroup::kspace_contribution() { double *vector_kspace = force->kspace->f2group; force->kspace->compute_group_group(groupbit,jgroupbit,0); scalar += 2.0*force->kspace->e2group; vector[0] += vector_kspace[0]; vector[1] += vector_kspace[1]; vector[2] += vector_kspace[2]; // subtract extra A <--> A Kspace interaction so energy matches // real-space style of compute group-group // add extra Kspace term to energy force->kspace->compute_group_group(groupbit,jgroupbit,1); scalar -= force->kspace->e2group; // self energy correction term scalar -= e_self; // k=0 boundary correction term if (boundaryflag) { double xprd = domain->xprd; double yprd = domain->yprd; double zprd = domain->zprd; // adjustment of z dimension for 2d slab Ewald // 3d Ewald just uses zprd since slab_volfactor = 1.0 double volume = xprd*yprd*zprd*force->kspace->slab_volfactor; scalar -= e_correction/volume; } } /* ---------------------------------------------------------------------- */ void ComputeGroupGroup::kspace_correction() { // total charge of groups A & B, needed for correction term double qsqsum_group,qsum_A,qsum_B; qsqsum_group = qsum_A = qsum_B = 0.0; double *q = atom->q; int *mask = atom->mask; int groupbit_A = groupbit; int groupbit_B = jgroupbit; for (int i = 0; i < atom->nlocal; i++) { if ((mask[i] & groupbit_A) && (mask[i] & groupbit_B)) qsqsum_group += q[i]*q[i]; if (mask[i] & groupbit_A) qsum_A += q[i]; if (mask[i] & groupbit_B) qsum_B += q[i]; } double tmp; MPI_Allreduce(&qsqsum_group,&tmp,1,MPI_DOUBLE,MPI_SUM,world); qsqsum_group = tmp; MPI_Allreduce(&qsum_A,&tmp,1,MPI_DOUBLE,MPI_SUM,world); qsum_A = tmp; MPI_Allreduce(&qsum_B,&tmp,1,MPI_DOUBLE,MPI_SUM,world); qsum_B = tmp; double g_ewald = force->kspace->g_ewald; double scale = 1.0; const double qscale = force->qqrd2e * scale; // self-energy correction e_self = qscale * g_ewald*qsqsum_group/MY_PIS; e_correction = 2.0*qsum_A*qsum_B; // subtract extra AA terms qsum_A = qsum_B = 0.0; for (int i = 0; i < atom->nlocal; i++) { if (!((mask[i] & groupbit_A) && (mask[i] & groupbit_B))) continue; if (mask[i] & groupbit_A) qsum_A += q[i]; if (mask[i] & groupbit_B) qsum_B += q[i]; } MPI_Allreduce(&qsum_A,&tmp,1,MPI_DOUBLE,MPI_SUM,world); qsum_A = tmp; MPI_Allreduce(&qsum_B,&tmp,1,MPI_DOUBLE,MPI_SUM,world); qsum_B = tmp; // k=0 energy correction term (still need to divide by volume above) e_correction -= qsum_A*qsum_B; e_correction *= qscale * MY_PI2 / (g_ewald*g_ewald); } diff --git a/src/compute_group_group.h b/src/compute_group_group.h index c881b08cb..906a185cf 100644 --- a/src/compute_group_group.h +++ b/src/compute_group_group.h @@ -1,89 +1,89 @@ /* -*- 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. ------------------------------------------------------------------------- */ #ifdef COMPUTE_CLASS ComputeStyle(group/group,ComputeGroupGroup) #else #ifndef LMP_COMPUTE_GROUP_GROUP_H #define LMP_COMPUTE_GROUP_GROUP_H #include "compute.h" namespace LAMMPS_NS { class ComputeGroupGroup : public Compute { public: ComputeGroupGroup(class LAMMPS *, int, char **); ~ComputeGroupGroup(); void init(); void init_list(int, class NeighList *); double compute_scalar(); void compute_vector(); private: char *group2; int jgroup,jgroupbit,othergroupbit; double **cutsq; double e_self,e_correction; - int pairflag,kspaceflag,boundaryflag; + int pairflag,kspaceflag,boundaryflag,molflag; class Pair *pair; class NeighList *list; class KSpace *kspace; void pair_contribution(); void kspace_contribution(); void kspace_correction(); }; } #endif #endif /* ERROR/WARNING messages: E: Illegal ... command Self-explanatory. Check the input script syntax and compare to the documentation for the command. You can use -echo screen as a command-line option when running LAMMPS to see the offending line. E: Compute group/group group ID does not exist Self-explanatory. E: No pair style defined for compute group/group Cannot calculate group interactions without a pair style defined. E: Pair style does not support compute group/group The pair_style does not have a single() function, so it cannot be invokded by the compute group/group command. E: No Kspace style defined for compute group/group Self-explanatory. E: Kspace style does not support compute group/group Self-explanatory. W: Both groups in compute group/group have a net charge; the Kspace boundary correction to energy will be non-zero Self-explantory. */ diff --git a/src/fix.cpp b/src/fix.cpp index 8c128b96b..0a95bcc69 100644 --- a/src/fix.cpp +++ b/src/fix.cpp @@ -1,310 +1,311 @@ /* ---------------------------------------------------------------------- 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. ------------------------------------------------------------------------- */ #include #include #include "fix.h" #include "atom.h" #include "group.h" #include "force.h" #include "comm.h" #include "atom_masks.h" #include "memory.h" #include "error.h" using namespace LAMMPS_NS; using namespace FixConst; // allocate space for static class instance variable and initialize it int Fix::instance_total = 0; /* ---------------------------------------------------------------------- */ -Fix::Fix(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp), -id(NULL), style(NULL), extlist(NULL), vector_atom(NULL), array_atom(NULL), -vector_local(NULL), array_local(NULL), eatom(NULL), vatom(NULL) +Fix::Fix(LAMMPS *lmp, int narg, char **arg) : + Pointers(lmp), + id(NULL), style(NULL), extlist(NULL), vector_atom(NULL), array_atom(NULL), + vector_local(NULL), array_local(NULL), eatom(NULL), vatom(NULL) { instance_me = instance_total++; // fix ID, group, and style // ID must be all alphanumeric chars or underscores int n = strlen(arg[0]) + 1; id = new char[n]; strcpy(id,arg[0]); for (int i = 0; i < n-1; i++) if (!isalnum(id[i]) && id[i] != '_') error->all(FLERR,"Fix ID must be alphanumeric or underscore characters"); igroup = group->find(arg[1]); if (igroup == -1) error->all(FLERR,"Could not find fix group ID"); groupbit = group->bitmask[igroup]; n = strlen(arg[2]) + 1; style = new char[n]; strcpy(style,arg[2]); restart_global = restart_peratom = restart_file = 0; force_reneighbor = 0; box_change_size = box_change_shape = box_change_domain = 0; thermo_energy = 0; rigid_flag = 0; peatom_flag = 0; virial_flag = 0; no_change_box = 0; time_integrate = 0; time_depend = 0; create_attribute = 0; restart_pbc = 0; wd_header = wd_section = 0; dynamic_group_allow = 0; dof_flag = 0; special_alter_flag = 0; enforce2d_flag = 0; respa_level_support = 0; respa_level = -1; scalar_flag = vector_flag = array_flag = 0; peratom_flag = local_flag = 0; size_vector_variable = size_array_rows_variable = 0; comm_forward = comm_reverse = comm_border = 0; restart_reset = 0; // reasonable defaults // however, each fix that uses these values should explicitly set them nevery = 1; global_freq = 1; // per-atom virial // set vflag_atom = 0 b/c some fixes grow vatom in grow_arrays() // which may occur outside of timestepping maxeatom = maxvatom = 0; vflag_atom = 0; // KOKKOS per-fix data masks execution_space = Host; datamask_read = ALL_MASK; datamask_modify = ALL_MASK; kokkosable = 0; copymode = 0; } /* ---------------------------------------------------------------------- */ Fix::~Fix() { if (copymode) return; delete [] id; delete [] style; memory->destroy(eatom); memory->destroy(vatom); } /* ---------------------------------------------------------------------- process params common to all fixes here if unknown param, call modify_param specific to the fix ------------------------------------------------------------------------- */ void Fix::modify_params(int narg, char **arg) { if (narg == 0) error->all(FLERR,"Illegal fix_modify command"); int iarg = 0; while (iarg < narg) { if (strcmp(arg[iarg],"energy") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix_modify command"); if (strcmp(arg[iarg+1],"no") == 0) thermo_energy = 0; else if (strcmp(arg[iarg+1],"yes") == 0) thermo_energy = 1; else error->all(FLERR,"Illegal fix_modify command"); iarg += 2; } else if (strcmp(arg[iarg],"respa") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix_modify command"); if (!respa_level_support) error->all(FLERR,"Illegal fix_modify command"); int lvl = force->inumeric(FLERR,arg[iarg+1]); if (lvl < 0) error->all(FLERR,"Illegal fix_modify command"); respa_level = lvl-1; iarg += 2; } else { int n = modify_param(narg-iarg,&arg[iarg]); if (n == 0) error->all(FLERR,"Illegal fix_modify command"); iarg += n; } } } /* ---------------------------------------------------------------------- setup for energy, virial computation see integrate::ev_set() for values of eflag (0-3) and vflag (0-6) fixes call this if use ev_tally() ------------------------------------------------------------------------- */ void Fix::ev_setup(int eflag, int vflag) { int i,n; evflag = 1; eflag_either = eflag; eflag_global = eflag % 2; eflag_atom = eflag / 2; vflag_either = vflag; vflag_global = vflag % 4; vflag_atom = vflag / 4; // reallocate per-atom arrays if necessary if (eflag_atom && atom->nlocal > maxeatom) { maxeatom = atom->nmax; memory->destroy(eatom); memory->create(eatom,maxeatom,"fix:eatom"); } if (vflag_atom && atom->nlocal > maxvatom) { maxvatom = atom->nmax; memory->destroy(vatom); memory->create(vatom,maxvatom,6,"fix:vatom"); } // zero accumulators // no global energy variable to zero (unlike pair,bond,angle,etc) // fixes tally it individually via fix_modify energy yes and compute_scalar() if (vflag_global) for (i = 0; i < 6; i++) virial[i] = 0.0; if (eflag_atom) { n = atom->nlocal; for (i = 0; i < n; i++) eatom[i] = 0.0; } if (vflag_atom) { n = atom->nlocal; for (i = 0; i < n; i++) { vatom[i][0] = 0.0; vatom[i][1] = 0.0; vatom[i][2] = 0.0; vatom[i][3] = 0.0; vatom[i][4] = 0.0; vatom[i][5] = 0.0; } } } /* ---------------------------------------------------------------------- setup for virial computation see integrate::ev_set() for values of vflag (0-6) fixes call this if use v_tally() ------------------------------------------------------------------------- */ void Fix::v_setup(int vflag) { int i,n; evflag = 1; vflag_global = vflag % 4; vflag_atom = vflag / 4; // reallocate per-atom array if necessary if (vflag_atom && atom->nlocal > maxvatom) { maxvatom = atom->nmax; memory->destroy(vatom); memory->create(vatom,maxvatom,6,"fix:vatom"); } // zero accumulators if (vflag_global) for (i = 0; i < 6; i++) virial[i] = 0.0; if (vflag_atom) { n = atom->nlocal; for (i = 0; i < n; i++) { vatom[i][0] = 0.0; vatom[i][1] = 0.0; vatom[i][2] = 0.0; vatom[i][3] = 0.0; vatom[i][4] = 0.0; vatom[i][5] = 0.0; } } } /* ---------------------------------------------------------------------- tally per-atom energy and global/per-atom virial into accumulators n = # of local owned atoms involved, with local indices in list eng = total energy for the interaction involving total atoms v = total virial for the interaction involving total atoms increment per-atom energy of each atom in list by 1/total fraction v_tally tallies virial this method can be used when fix computes energy/forces in post_force() e.g. fix cmap: compute energy and virial only on owned atoms whether newton_bond is on or off other procs will tally left-over fractions for atoms they own ------------------------------------------------------------------------- */ void Fix::ev_tally(int n, int *list, double total, double eng, double *v) { if (eflag_atom) { double fraction = eng/total; for (int i = 0; i < n; i++) eatom[list[i]] += fraction; } v_tally(n,list,total,v); } /* ---------------------------------------------------------------------- tally virial into global and per-atom accumulators n = # of local owned atoms involved, with local indices in list v = total virial for the interaction involving total atoms increment global virial by n/total fraction increment per-atom virial of each atom in list by 1/total fraction this method can be used when fix computes forces in post_force() e.g. fix shake, fix rigid: compute virial only on owned atoms whether newton_bond is on or off other procs will tally left-over fractions for atoms they own ------------------------------------------------------------------------- */ void Fix::v_tally(int n, int *list, double total, double *v) { int m; if (vflag_global) { double fraction = n/total; virial[0] += fraction*v[0]; virial[1] += fraction*v[1]; virial[2] += fraction*v[2]; virial[3] += fraction*v[3]; virial[4] += fraction*v[4]; virial[5] += fraction*v[5]; } if (vflag_atom) { double fraction = 1.0/total; for (int i = 0; i < n; i++) { m = list[i]; vatom[m][0] += fraction*v[0]; vatom[m][1] += fraction*v[1]; vatom[m][2] += fraction*v[2]; vatom[m][3] += fraction*v[3]; vatom[m][4] += fraction*v[4]; vatom[m][5] += fraction*v[5]; } } } diff --git a/src/pair_buck.cpp b/src/pair_buck.cpp index 78a5321cf..ac15e8202 100644 --- a/src/pair_buck.cpp +++ b/src/pair_buck.cpp @@ -1,408 +1,409 @@ /* ---------------------------------------------------------------------- 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. ------------------------------------------------------------------------- */ #include #include #include #include #include "pair_buck.h" #include "atom.h" #include "comm.h" #include "force.h" #include "neigh_list.h" #include "math_const.h" #include "memory.h" #include "error.h" using namespace LAMMPS_NS; using namespace MathConst; /* ---------------------------------------------------------------------- */ PairBuck::PairBuck(LAMMPS *lmp) : Pair(lmp) { writedata = 1; } /* ---------------------------------------------------------------------- */ PairBuck::~PairBuck() { if (copymode) return; if (allocated) { memory->destroy(setflag); memory->destroy(cutsq); memory->destroy(cut); memory->destroy(a); memory->destroy(rho); memory->destroy(c); memory->destroy(rhoinv); memory->destroy(buck1); memory->destroy(buck2); memory->destroy(offset); } } /* ---------------------------------------------------------------------- */ void PairBuck::compute(int eflag, int vflag) { int i,j,ii,jj,inum,jnum,itype,jtype; double xtmp,ytmp,ztmp,delx,dely,delz,evdwl,fpair; double rsq,r2inv,r6inv,forcebuck,factor_lj; double r,rexp; int *ilist,*jlist,*numneigh,**firstneigh; evdwl = 0.0; if (eflag || vflag) ev_setup(eflag,vflag); else evflag = vflag_fdotr = 0; double **x = atom->x; double **f = atom->f; int *type = atom->type; int nlocal = atom->nlocal; double *special_lj = force->special_lj; int newton_pair = force->newton_pair; inum = list->inum; ilist = list->ilist; numneigh = list->numneigh; firstneigh = list->firstneigh; // loop over neighbors of my atoms for (ii = 0; ii < inum; ii++) { i = ilist[ii]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; itype = type[i]; jlist = firstneigh[i]; jnum = numneigh[i]; for (jj = 0; jj < jnum; jj++) { j = jlist[jj]; factor_lj = special_lj[sbmask(j)]; j &= NEIGHMASK; delx = xtmp - x[j][0]; dely = ytmp - x[j][1]; delz = ztmp - x[j][2]; rsq = delx*delx + dely*dely + delz*delz; jtype = type[j]; if (rsq < cutsq[itype][jtype]) { r2inv = 1.0/rsq; r6inv = r2inv*r2inv*r2inv; r = sqrt(rsq); rexp = exp(-r*rhoinv[itype][jtype]); forcebuck = buck1[itype][jtype]*r*rexp - buck2[itype][jtype]*r6inv; fpair = factor_lj*forcebuck*r2inv; f[i][0] += delx*fpair; f[i][1] += dely*fpair; f[i][2] += delz*fpair; if (newton_pair || j < nlocal) { f[j][0] -= delx*fpair; f[j][1] -= dely*fpair; f[j][2] -= delz*fpair; } if (eflag) { evdwl = a[itype][jtype]*rexp - c[itype][jtype]*r6inv - offset[itype][jtype]; evdwl *= factor_lj; } if (evflag) ev_tally(i,j,nlocal,newton_pair, evdwl,0.0,fpair,delx,dely,delz); } } } if (vflag_fdotr) virial_fdotr_compute(); } /* ---------------------------------------------------------------------- allocate all arrays ------------------------------------------------------------------------- */ void PairBuck::allocate() { allocated = 1; int n = atom->ntypes; memory->create(setflag,n+1,n+1,"pair:setflag"); for (int i = 1; i <= n; i++) for (int j = i; j <= n; j++) setflag[i][j] = 0; memory->create(cutsq,n+1,n+1,"pair:cutsq"); memory->create(cut,n+1,n+1,"pair:cut_lj"); memory->create(a,n+1,n+1,"pair:a"); memory->create(rho,n+1,n+1,"pair:rho"); memory->create(c,n+1,n+1,"pair:c"); memory->create(rhoinv,n+1,n+1,"pair:rhoinv"); memory->create(buck1,n+1,n+1,"pair:buck1"); memory->create(buck2,n+1,n+1,"pair:buck2"); memory->create(offset,n+1,n+1,"pair:offset"); } /* ---------------------------------------------------------------------- global settings ------------------------------------------------------------------------- */ void PairBuck::settings(int narg, char **arg) { if (narg != 1) error->all(FLERR,"Illegal pair_style command"); cut_global = force->numeric(FLERR,arg[0]); // reset cutoffs that have been explicitly set if (allocated) { int i,j; for (i = 1; i <= atom->ntypes; i++) for (j = i+1; j <= atom->ntypes; j++) if (setflag[i][j]) cut[i][j] = cut_global; } } /* ---------------------------------------------------------------------- set coeffs for one or more type pairs ------------------------------------------------------------------------- */ void PairBuck::coeff(int narg, char **arg) { - if (narg < 5 || narg > 6) error->all(FLERR,"Incorrect args for pair coefficients"); + if (narg < 5 || narg > 6) + error->all(FLERR,"Incorrect args for pair coefficients"); if (!allocated) allocate(); int ilo,ihi,jlo,jhi; force->bounds(FLERR,arg[0],atom->ntypes,ilo,ihi); force->bounds(FLERR,arg[1],atom->ntypes,jlo,jhi); double a_one = force->numeric(FLERR,arg[2]); double rho_one = force->numeric(FLERR,arg[3]); if (rho_one <= 0) error->all(FLERR,"Incorrect args for pair coefficients"); double c_one = force->numeric(FLERR,arg[4]); double cut_one = cut_global; if (narg == 6) cut_one = force->numeric(FLERR,arg[5]); int count = 0; for (int i = ilo; i <= ihi; i++) { for (int j = MAX(jlo,i); j <= jhi; j++) { a[i][j] = a_one; rho[i][j] = rho_one; c[i][j] = c_one; cut[i][j] = cut_one; setflag[i][j] = 1; count++; } } if (count == 0) error->all(FLERR,"Incorrect args for pair coefficients"); } /* ---------------------------------------------------------------------- init for one type pair i,j and corresponding j,i ------------------------------------------------------------------------- */ double PairBuck::init_one(int i, int j) { if (setflag[i][j] == 0) error->all(FLERR,"All pair coeffs are not set"); rhoinv[i][j] = 1.0/rho[i][j]; buck1[i][j] = a[i][j]/rho[i][j]; buck2[i][j] = 6.0*c[i][j]; if (offset_flag) { double rexp = exp(-cut[i][j]/rho[i][j]); offset[i][j] = a[i][j]*rexp - c[i][j]/pow(cut[i][j],6.0); } else offset[i][j] = 0.0; a[j][i] = a[i][j]; c[j][i] = c[i][j]; rhoinv[j][i] = rhoinv[i][j]; buck1[j][i] = buck1[i][j]; buck2[j][i] = buck2[i][j]; offset[j][i] = offset[i][j]; // compute I,J contribution to long-range tail correction // count total # of atoms of type I and J via Allreduce if (tail_flag) { int *type = atom->type; int nlocal = atom->nlocal; double count[2],all[2]; count[0] = count[1] = 0.0; for (int k = 0; k < nlocal; k++) { if (type[k] == i) count[0] += 1.0; if (type[k] == j) count[1] += 1.0; } MPI_Allreduce(count,all,2,MPI_DOUBLE,MPI_SUM,world); double rho1 = rho[i][j]; double rho2 = rho1*rho1; double rho3 = rho2*rho1; double rc = cut[i][j]; double rc2 = rc*rc; double rc3 = rc2*rc; etail_ij = 2.0*MY_PI*all[0]*all[1]* (a[i][j]*exp(-rc/rho1)*rho1*(rc2 + 2.0*rho1*rc + 2.0*rho2) - c[i][j]/(3.0*rc3)); ptail_ij = (-1/3.0)*2.0*MY_PI*all[0]*all[1]* (-a[i][j]*exp(-rc/rho1)* (rc3 + 3.0*rho1*rc2 + 6.0*rho2*rc + 6.0*rho3) + 2.0*c[i][j]/rc3); } return cut[i][j]; } /* ---------------------------------------------------------------------- proc 0 writes to restart file ------------------------------------------------------------------------- */ void PairBuck::write_restart(FILE *fp) { write_restart_settings(fp); int i,j; for (i = 1; i <= atom->ntypes; i++) for (j = i; j <= atom->ntypes; j++) { fwrite(&setflag[i][j],sizeof(int),1,fp); if (setflag[i][j]) { fwrite(&a[i][j],sizeof(double),1,fp); fwrite(&rho[i][j],sizeof(double),1,fp); fwrite(&c[i][j],sizeof(double),1,fp); fwrite(&cut[i][j],sizeof(double),1,fp); } } } /* ---------------------------------------------------------------------- proc 0 reads from restart file, bcasts ------------------------------------------------------------------------- */ void PairBuck::read_restart(FILE *fp) { read_restart_settings(fp); allocate(); int i,j; int me = comm->me; for (i = 1; i <= atom->ntypes; i++) for (j = i; j <= atom->ntypes; j++) { if (me == 0) fread(&setflag[i][j],sizeof(int),1,fp); MPI_Bcast(&setflag[i][j],1,MPI_INT,0,world); if (setflag[i][j]) { if (me == 0) { fread(&a[i][j],sizeof(double),1,fp); fread(&rho[i][j],sizeof(double),1,fp); fread(&c[i][j],sizeof(double),1,fp); fread(&cut[i][j],sizeof(double),1,fp); } MPI_Bcast(&a[i][j],1,MPI_DOUBLE,0,world); MPI_Bcast(&rho[i][j],1,MPI_DOUBLE,0,world); MPI_Bcast(&c[i][j],1,MPI_DOUBLE,0,world); MPI_Bcast(&cut[i][j],1,MPI_DOUBLE,0,world); } } } /* ---------------------------------------------------------------------- proc 0 writes to restart file ------------------------------------------------------------------------- */ void PairBuck::write_restart_settings(FILE *fp) { fwrite(&cut_global,sizeof(double),1,fp); fwrite(&offset_flag,sizeof(int),1,fp); fwrite(&mix_flag,sizeof(int),1,fp); fwrite(&tail_flag,sizeof(int),1,fp); } /* ---------------------------------------------------------------------- proc 0 reads from restart file, bcasts ------------------------------------------------------------------------- */ void PairBuck::read_restart_settings(FILE *fp) { if (comm->me == 0) { fread(&cut_global,sizeof(double),1,fp); fread(&offset_flag,sizeof(int),1,fp); fread(&mix_flag,sizeof(int),1,fp); fread(&tail_flag,sizeof(int),1,fp); } MPI_Bcast(&cut_global,1,MPI_DOUBLE,0,world); MPI_Bcast(&offset_flag,1,MPI_INT,0,world); MPI_Bcast(&mix_flag,1,MPI_INT,0,world); MPI_Bcast(&tail_flag,1,MPI_INT,0,world); } /* ---------------------------------------------------------------------- proc 0 writes to data file ------------------------------------------------------------------------- */ void PairBuck::write_data(FILE *fp) { for (int i = 1; i <= atom->ntypes; i++) fprintf(fp,"%d %g %g %g\n",i,a[i][i],rho[i][i],c[i][i]); } /* ---------------------------------------------------------------------- proc 0 writes all pairs to data file ------------------------------------------------------------------------- */ void PairBuck::write_data_all(FILE *fp) { for (int i = 1; i <= atom->ntypes; i++) for (int j = i; j <= atom->ntypes; j++) fprintf(fp,"%d %d %g %g %g %g\n",i,j, a[i][j],rho[i][j],c[i][j],cut[i][j]); } /* ---------------------------------------------------------------------- */ double PairBuck::single(int i, int j, int itype, int jtype, double rsq, double factor_coul, double factor_lj, double &fforce) { double r2inv,r6inv,r,rexp,forcebuck,phibuck; r2inv = 1.0/rsq; r6inv = r2inv*r2inv*r2inv; r = sqrt(rsq); rexp = exp(-r*rhoinv[itype][jtype]); forcebuck = buck1[itype][jtype]*r*rexp - buck2[itype][jtype]*r6inv; fforce = factor_lj*forcebuck*r2inv; phibuck = a[itype][jtype]*rexp - c[itype][jtype]*r6inv - offset[itype][jtype]; return factor_lj*phibuck; } /* ---------------------------------------------------------------------- */ void *PairBuck::extract(const char *str, int &dim) { dim = 2; if (strcmp(str,"a") == 0) return (void *) a; if (strcmp(str,"c") == 0) return (void *) c; return NULL; } diff --git a/src/version.h b/src/version.h index 6847e5ff8..299704ccd 100644 --- a/src/version.h +++ b/src/version.h @@ -1 +1 @@ -#define LAMMPS_VERSION "26 Jan 2017" +#define LAMMPS_VERSION "10 Feb 2016"