diff --git a/doc/src/Section_commands.txt b/doc/src/Section_commands.txt index 78969bc4a..32b714e1d 100644 --- a/doc/src/Section_commands.txt +++ b/doc/src/Section_commands.txt @@ -1,1247 +1,1248 @@ "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 core LAMMPS commands, grouped by category. The "next section"_#cmd_5 lists all 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 the categories here, 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 diagnostics (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 netcdf"_dump_netcdf.html, "dump netcdf/mpiio"_dump_netcdf.html, "dump vtk"_dump_vtk.html, "group2ndx"_group2ndx.html, "ndx2group"_group2ndx.html, "temper/grem"_temper_grem.html "temper/npt"_temper_npt.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, "python"_fix_python.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, "edpd/source"_fix_dpd_source.html, "eos/cv"_fix_eos_cv.html, "eos/table"_fix_eos_table.html, "eos/table/rx"_fix_eos_table_rx.html, "filter/corotate"_fix_filter_corotate.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, "mvv/dpd"_fix_mvv_dpd.html, "mvv/edpd"_fix_mvv_dpd.html, "mvv/tdpd"_fix_mvv_dpd.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 (ko)"_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, "tdpd/source"_fix_dpd_source.html, "temp/rescale/eff"_fix_temp_rescale_eff.html, "ti/spring"_fix_ti_spring.html, "ttm/mod"_fix_ttm.html, "wall/ees"_fix_wall_ees.html, "wall/region/ees"_fix_wall_ees.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. +"aggregate/atom"_compute_cluster_atom.html, "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, "fragment/atom"_compute_cluster_atom.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, "cnp/atom"_compute_cnp_atom.html, "dpd"_compute_dpd.html, "dpd/atom"_compute_dpd_atom.html, "edpd/temp/atom"_compute_edpd_temp_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, "tdpd/cc/atom"_compute_tdpd_cc_atom.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 (oi)"_pair_airebo.html, "airebo/morse (oi)"_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 (gkiot)"_pair_eam.html, "eam/fs (gkiot)"_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, "gw"_pair_gw.html, "gw/zbl"_pair_gw.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 (kio)"_pair_charmm.html, "lj/charmm/coul/charmm/implicit (ko)"_pair_charmm.html, "lj/charmm/coul/long (gkio)"_pair_charmm.html, "lj/charmm/coul/msm"_pair_charmm.html, "lj/charmmfsw/coul/charmmfsh"_pair_charmm.html, "lj/charmmfsw/coul/long"_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 (io)"_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 (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, "python"_pair_python.html, "reax"_pair_reax.html, "rebo (oi)"_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, "edip/multi"_pair_edip.html, "edpd"_pair_meso.html, "eff/cut"_pair_eff.html, "exp6/rx"_pair_exp6_rx.html, "gauss/cut"_pair_gauss.html, "kolmogorov/crespi/z"_pair_kolmogorov_crespi_z.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, "mdpd"_pair_meso.html, "mdpd/rhosum"_pair_meso.html, "meam/c"_pair_meam.html, "meam/spline (o)"_pair_meam_spline.html, "meam/sw/spline"_pair_meam_sw_spline.html, "mgpt"_pair_mgpt.html, "momb"_pair_momb.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, "oxdna2/coaxstk"_pair_oxdna2.html, "oxdna2/dh"_pair_oxdna2.html, "oxdna2/excv"_pair_oxdna2.html, "oxdna2/stk"_pair_oxdna2.html, "quip"_pair_quip.html, "reax/c (ko)"_pair_reaxc.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, "tdpd"_pair_meso.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 (ko)"_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.html, "oxdna2/fene"_bond_oxdna.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 (ko)"_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, "charmmfsw"_dihedral_charmm.html, "class2 (ko)"_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 (ko)"_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 (i)"_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.txt b/doc/src/compute.txt index eed190aeb..c06735d28 100644 --- a/doc/src/compute.txt +++ b/doc/src/compute.txt @@ -1,262 +1,264 @@ "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 command :h3 [Syntax:] compute ID group-ID style args :pre ID = user-assigned name for the computation group-ID = ID of the group of atoms to perform the computation on style = one of a list of possible style names (see below) args = arguments used by a particular style :ul [Examples:] compute 1 all temp compute newtemp flow temp/partial 1 1 0 compute 3 all ke/atom :pre [Description:] Define a computation that will be performed on a group of atoms. Quantities calculated by a compute are instantaneous values, meaning they are calculated from information about atoms on the current timestep or iteration, though a compute may internally store some information about a previous state of the system. Defining a compute does not perform a computation. Instead computes are invoked by other LAMMPS commands as needed, e.g. to calculate a temperature needed for a thermostat fix or to generate thermodynamic or dump file output. See this "howto section"_Section_howto.html#howto_15 for a summary of various LAMMPS output options, many of which involve computes. The ID of a compute can only contain alphanumeric characters and underscores. :line Computes calculate one of three styles of quantities: global, per-atom, or local. A global quantity is one or more system-wide values, e.g. the temperature of the system. A per-atom quantity is one or more values per atom, e.g. the kinetic energy of each atom. Per-atom values are set to 0.0 for atoms not in the specified compute group. Local quantities are calculated by each processor based on the atoms it owns, but there may be zero or more per atom, e.g. a list of bond distances. Computes that produce per-atom quantities have the word "atom" in their style, e.g. {ke/atom}. Computes that produce local quantities have the word "local" in their style, e.g. {bond/local}. Styles with neither "atom" or "local" in their style produce global quantities. Note that a single compute produces either global or per-atom or local quantities, but never more than one of these (with only a few exceptions, as documented by individual compute commands). Global, per-atom, and local quantities each come in three kinds: a single scalar value, a vector of values, or a 2d array of values. The doc page for each compute describes the style and kind of values it produces, e.g. a per-atom vector. Some computes produce more than one kind of a single style, e.g. a global scalar and a global vector. When a compute quantity is accessed, as in many of the output commands discussed below, it can be referenced via the following bracket notation, where ID is the ID of the compute: c_ID | entire scalar, vector, or array c_ID\[I\] | one element of vector, one column of array c_ID\[I\]\[J\] | one element of array :tb(s=|) In other words, using one bracket reduces the dimension of the quantity once (vector -> scalar, array -> vector). Using two brackets reduces the dimension twice (array -> scalar). Thus a command that uses scalar compute values as input can also process elements of a vector or array. Note that commands and "variables"_variable.html which use compute quantities typically do not allow for all kinds, e.g. a command may require a vector of values, not a scalar. This means there is no ambiguity about referring to a compute quantity as c_ID even if it produces, for example, both a scalar and vector. The doc pages for various commands explain the details. :line In LAMMPS, the values generated by a compute can be used in several ways: The results of computes that calculate a global temperature or pressure can be used by fixes that do thermostatting or barostatting or when atom velocities are created. :ulb,l Global values can be output via the "thermo_style custom"_thermo_style.html or "fix ave/time"_fix_ave_time.html command. Or the values can be referenced in a "variable equal"_variable.html or "variable atom"_variable.html command. :l Per-atom values can be output via the "dump custom"_dump.html command. Or they can be time-averaged via the "fix ave/atom"_fix_ave_atom.html command or reduced by the "compute reduce"_compute_reduce.html command. Or the per-atom values can be referenced in an "atom-style variable"_variable.html. :l Local values can be reduced by the "compute reduce"_compute_reduce.html command, or histogrammed by the "fix ave/histo"_fix_ave_histo.html command, or output by the "dump local"_dump.html command. :l :ule The results of computes that calculate global quantities can be either "intensive" or "extensive" values. Intensive means the value is independent of the number of atoms in the simulation, e.g. temperature. Extensive means the value scales with the number of atoms in the simulation, e.g. total rotational kinetic energy. "Thermodynamic output"_thermo_style.html will normalize extensive values by the number of atoms in the system, depending on the "thermo_modify norm" setting. It will not normalize intensive values. If a compute value is accessed in another way, e.g. by a "variable"_variable.html, you may want to know whether it is an intensive or extensive value. See the doc page for individual computes for further info. :line LAMMPS creates its own computes internally for thermodynamic output. Three computes are always created, named "thermo_temp", "thermo_press", and "thermo_pe", as if these commands had been invoked in the input script: compute thermo_temp all temp compute thermo_press all pressure thermo_temp compute thermo_pe all pe :pre Additional computes for other quantities are created if the thermo style requires it. See the documentation for the "thermo_style"_thermo_style.html command. Fixes that calculate temperature or pressure, i.e. for thermostatting or barostatting, may also create computes. These are discussed in the documentation for specific "fix"_fix.html commands. In all these cases, the default computes LAMMPS creates can be replaced by computes defined by the user in the input script, as described by the "thermo_modify"_thermo_modify.html and "fix modify"_fix_modify.html commands. Properties of either a default or user-defined compute can be modified via the "compute_modify"_compute_modify.html command. Computes can be deleted with the "uncompute"_uncompute.html command. Code for new computes can be added to LAMMPS (see "this section"_Section_modify.html of the manual) and the results of their calculations accessed in the various ways described above. :line Each compute style has its own doc page which describes its arguments and what it does. Here is an alphabetic list of compute styles available in LAMMPS. They are also given in more compact form in the Compute section of "this page"_Section_commands.html#cmd_5. There are also additional compute styles (not listed here) submitted by users which are included in the LAMMPS distribution. The list of these with links to the individual styles are given in the compute section of "this page"_Section_commands.html#cmd_5. +"aggregate/atom"_compute_cluster_atom.html - aggregate ID for each atom "angle/local"_compute_bond_local.html - theta and energy of each angle "angmom/chunk"_compute_angmom_chunk.html - angular momentum for each chunk "body/local"_compute_body_local.html - attributes of body sub-particles "bond"_compute_bond.html - values computed by a bond style "bond/local"_compute_bond_local.html - distance and energy of each bond "centro/atom"_compute_centro_atom.html - centro-symmetry parameter for each atom "chunk/atom"_compute_chunk_atom.html - assign chunk IDs to each atom "cluster/atom"_compute_cluster_atom.html - cluster ID for each atom "cna/atom"_compute_cna_atom.html - common neighbor analysis (CNA) for each atom "com"_compute_com.html - center-of-mass of group of atoms "com/chunk"_compute_com_chunk.html - center-of-mass for each chunk "contact/atom"_compute_contact_atom.html - contact count for each spherical particle "coord/atom"_compute_coord_atom.html - coordination number for each atom "damage/atom"_compute_damage_atom.html - Peridynamic damage for each atom "dihedral/local"_compute_dihedral_local.html - angle of each dihedral "dilatation/atom"_compute_dilatation_atom.html - Peridynamic dilatation for each atom "displace/atom"_compute_displace_atom.html - displacement of each atom "erotate/asphere"_compute_erotate_asphere.html - rotational energy of aspherical particles "erotate/rigid"_compute_erotate_rigid.html - rotational energy of rigid bodies "erotate/sphere"_compute_erotate_sphere.html - rotational energy of spherical particles "erotate/sphere/atom"_compute_erotate_sphere.html - rotational energy for each spherical particle "event/displace"_compute_event_displace.html - detect event on atom displacement +"fragment/atom"_compute_cluster_atom.html - fragment ID for each atom "group/group"_compute_group_group.html - energy/force between two groups of atoms "gyration"_compute_gyration.html - radius of gyration of group of atoms "gyration/chunk"_compute_gyration_chunk.html - radius of gyration for each chunk "heat/flux"_compute_heat_flux.html - heat flux through a group of atoms "hexorder/atom"_compute_hexorder_atom.html - bond orientational order parameter q6 "improper/local"_compute_improper_local.html - angle of each improper "inertia/chunk"_compute_inertia_chunk.html - inertia tensor for each chunk "ke"_compute_ke.html - translational kinetic energy "ke/atom"_compute_ke_atom.html - kinetic energy for each atom "ke/rigid"_compute_ke_rigid.html - translational kinetic energy of rigid bodies "msd"_compute_msd.html - mean-squared displacement of group of atoms "msd/chunk"_compute_msd_chunk.html - mean-squared displacement for each chunk "msd/nongauss"_compute_msd_nongauss.html - MSD and non-Gaussian parameter of group of atoms "omega/chunk"_compute_omega_chunk.html - angular velocity for each chunk "orientorder/atom"_compute_orientorder_atom.html - Steinhardt bond orientational order parameters Ql "pair"_compute_pair.html - values computed by a pair style "pair/local"_compute_pair_local.html - distance/energy/force of each pairwise interaction "pe"_compute_pe.html - potential energy "pe/atom"_compute_pe_atom.html - potential energy for each atom "plasticity/atom"_compute_plasticity_atom.html - Peridynamic plasticity for each atom "pressure"_compute_pressure.html - total pressure and pressure tensor "property/atom"_compute_property_atom.html - convert atom attributes to per-atom vectors/arrays "property/local"_compute_property_local.html - convert local attributes to localvectors/arrays "property/chunk"_compute_property_chunk.html - extract various per-chunk attributes "rdf"_compute_rdf.html - radial distribution function g(r) histogram of group of atoms "reduce"_compute_reduce.html - combine per-atom quantities into a single global value "reduce/region"_compute_reduce.html - same as compute reduce, within a region "rigid/local"_compute_rigid_local.html - extract rigid body attributes "slice"_compute_slice.html - extract values from global vector or array "sna/atom"_compute_sna_atom.html - calculate bispectrum coefficients for each atom "snad/atom"_compute_sna_atom.html - derivative of bispectrum coefficients for each atom "snav/atom"_compute_sna_atom.html - virial contribution from bispectrum coefficients for each atom "stress/atom"_compute_stress_atom.html - stress tensor for each atom "temp"_compute_temp.html - temperature of group of atoms "temp/asphere"_compute_temp_asphere.html - temperature of aspherical particles "temp/body"_compute_temp_body.html - temperature of body particles "temp/chunk"_compute_temp_chunk.html - temperature of each chunk "temp/com"_compute_temp_com.html - temperature after subtracting center-of-mass velocity "temp/deform"_compute_temp_deform.html - temperature excluding box deformation velocity "temp/partial"_compute_temp_partial.html - temperature excluding one or more dimensions of velocity "temp/profile"_compute_temp_profile.html - temperature excluding a binned velocity profile "temp/ramp"_compute_temp_ramp.html - temperature excluding ramped velocity component "temp/region"_compute_temp_region.html - temperature of a region of atoms "temp/sphere"_compute_temp_sphere.html - temperature of spherical particles "ti"_compute_ti.html - thermodynamic integration free energy values "torque/chunk"_compute_torque_chunk.html - torque applied on each chunk "vacf"_compute_vacf.html - velocity-autocorrelation function of group of atoms "vcm/chunk"_compute_vcm_chunk.html - velocity of center-of-mass for each chunk "voronoi/atom"_compute_voronoi_atom.html - Voronoi volume and neighbors for each atom :ul There are also additional compute styles submitted by users which are included in the LAMMPS distribution. The list of these with links to the individual styles are given in the compute section of "this page"_Section_commands.html#cmd_5. There are also additional accelerated compute styles included in the LAMMPS distribution for faster performance on CPUs and GPUs. The list of these with links to the individual styles are given in the pair section of "this page"_Section_commands.html#cmd_5. [Restrictions:] none [Related commands:] "uncompute"_uncompute.html, "compute_modify"_compute_modify.html, "fix ave/atom"_fix_ave_atom.html, "fix ave/time"_fix_ave_time.html, "fix ave/histo"_fix_ave_histo.html [Default:] none diff --git a/doc/src/fix_saed_vtk.txt b/doc/src/fix_saed_vtk.txt index 3523233ce..814e17b8e 100644 --- a/doc/src/fix_saed_vtk.txt +++ b/doc/src/fix_saed_vtk.txt @@ -1,190 +1,190 @@ "LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c :link(lws,http://lammps.sandia.gov) :link(ld,Manual.html) :link(lc,Section_commands.html#comm) :line fix saed/vtk command :h3 [Syntax:] fix ID group-ID saed/vtk Nevery Nrepeat Nfreak c_ID attribute args ... keyword args ... :pre ID, group-ID are documented in "fix"_fix.html command :ulb,l saed/vtk = style name of this fix command :l Nevery = use input values every this many timesteps :l Nrepeat = # of times to use input values for calculating averages :l Nfreq = calculate averages every this many timesteps :l c_ID = saed compute ID :l keyword = {file} or {ave} or {start} or {file} or {overwrite}:l {ave} args = {one} or {running} or {window M} one = output a new average value every Nfreq steps running = output cumulative average of all previous Nfreq steps window M = output average of M most recent Nfreq steps {start} args = Nstart Nstart = start averaging on this timestep {file} arg = filename filename = name of file to output time averages to {overwrite} arg = none = overwrite output file with only latest output :pre :ule [Examples:] compute 1 all saed 0.0251 Al O Kmax 1.70 Zone 0 0 1 dR_Ewald 0.01 c 0.5 0.5 0.5 compute 2 all saed 0.0251 Ni Kmax 1.70 Zone 0 0 0 c 0.05 0.05 0.05 manual echo :pre -fix saed/vtk 1 1 1 c_1 file Al2O3_001.saed -fix saed/vtk 1 1 1 c_2 file Ni_000.saed :pre +fix 1 all saed/vtk 1 1 1 c_1 file Al2O3_001.saed +fix 2 all saed/vtk 1 1 1 c_2 file Ni_000.saed :pre [Description:] Time average computed intensities from "compute saed"_compute_saed.html and write output to a file in the 3rd generation vtk image data format for visualization directly in parallelized visualization software packages like ParaView and VisIt. Note that if no time averaging is done, this command can be used as a convenient way to simply output diffraction intensities at a single snapshot. To produce output in the image data vtk format ghost data is added outside the {Kmax} range assigned in the compute saed. The ghost data is assigned a value of -1 and can be removed setting a minimum isovolume of 0 within the visualization software. SAED images can be created by visualizing a spherical slice of the data that is centered at R_Ewald*\[h k l\]/norm(\[h k l\]), where R_Ewald=1/lambda. The group specified within this command is ignored. However, note that specified values may represent calculations performed by saed computes which store their own "group" definitions. Fix saed/vtk is designed to work only with "compute saed"_compute_saed.html values, e.g. compute 3 top saed 0.0251 Al O fix saed/vtk 1 1 1 c_3 file Al2O3_001.saed :pre :line The {Nevery}, {Nrepeat}, and {Nfreq} arguments specify on what timesteps the input values will be used in order to contribute to the average. The final averaged quantities are generated on timesteps that are a multiple of {Nfreq}. The average is over {Nrepeat} quantities, computed in the preceding portion of the simulation every {Nevery} timesteps. {Nfreq} must be a multiple of {Nevery} and {Nevery} must be non-zero even if {Nrepeat} is 1. Also, the timesteps contributing to the average value cannot overlap, i.e. Nrepeat*Nevery can not exceed Nfreq. For example, if Nevery=2, Nrepeat=6, and Nfreq=100, then values on timesteps 90,92,94,96,98,100 will be used to compute the final average on timestep 100. Similarly for timesteps 190,192,194,196,198,200 on timestep 200, etc. If Nrepeat=1 and Nfreq = 100, then no time averaging is done; values are simply generated on timesteps 100,200,etc. :line The output for fix ave/time/saed is a file written with the 3rd generation vtk image data formatting. The filename assigned by the {file} keyword is appended with _N.vtk where N is an index (0,1,2...) to account for multiple diffraction intensity outputs. By default the header contains the following information (with example data): # vtk DataFile Version 3.0 c_SAED Image data set ASCII DATASET STRUCTURED_POINTS DIMENSIONS 337 219 209 ASPECT_RATIO 0.00507953 0.00785161 0.00821458 ORIGIN -0.853361 -0.855826 -0.854316 POINT_DATA 15424827 SCALARS intensity float LOOKUP_TABLE default ...data :pre In this example, kspace is sampled across a 337 x 219 x 209 point mesh where the mesh spacing is approximately 0.005, 0.007, and 0.008 inv(length) units in the k1, k2, and k3 directions, respectively. The data is shifted by -0.85, -0.85, -0.85 inv(length) units so that the origin will lie at 0, 0, 0. Here, 15,424,827 kspace points are sampled in total. :line Additional optional keywords also affect the operation of this fix. The {ave} keyword determines how the values produced every {Nfreq} steps are averaged with values produced on previous steps that were multiples of {Nfreq}, before they are accessed by another output command or written to a file. If the {ave} setting is {one}, then the values produced on timesteps that are multiples of {Nfreq} are independent of each other; they are output as-is without further averaging. If the {ave} setting is {running}, then the values produced on timesteps that are multiples of {Nfreq} are summed and averaged in a cumulative sense before being output. Each output value is thus the average of the value produced on that timestep with all preceding values. This running average begins when the fix is defined; it can only be restarted by deleting the fix via the "unfix"_unfix.html command, or by re-defining the fix by re-specifying it. If the {ave} setting is {window}, then the values produced on timesteps that are multiples of {Nfreq} are summed and averaged within a moving "window" of time, so that the last M values are used to produce the output. E.g. if M = 3 and Nfreq = 1000, then the output on step 10000 will be the average of the individual values on steps 8000,9000,10000. Outputs on early steps will average over less than M values if they are not available. The {start} keyword specifies what timestep averaging will begin on. The default is step 0. Often input values can be 0.0 at time 0, so setting {start} to a larger value can avoid including a 0.0 in a running or windowed average. The {file} keyword allows a filename to be specified. Every {Nfreq} steps, the vector of saed intensity data is written to a new file using the 3rd generation vtk format. The base of each file is assigned by the {file} keyword and this string is appended with _N.vtk where N is an index (0,1,2...) to account for situations with multiple diffraction intensity outputs. The {overwrite} keyword will continuously overwrite the output file with the latest output, so that it only contains one timestep worth of output. This option can only be used with the {ave running} setting. [Restart, fix_modify, output, run start/stop, minimize info:] No information about this fix is written to "binary restart files"_restart.html. None of the "fix_modify"_fix_modify.html options are relevant to this fix. No parameter of this fix can be used with the {start/stop} keywords of the "run"_run.html command. This fix is not invoked during "energy minimization"_minimize.html. [Restrictions:] The attributes for fix_saed_vtk must match the values assigned in the associated "compute_saed"_compute_saed.html command. [Related commands:] "compute_saed"_compute_saed.html [Default:] The option defaults are ave = one, start = 0, no file output. :line :link(Coleman) [(Coleman)] Coleman, Spearot, Capolungo, MSMSE, 21, 055020 (2013). diff --git a/src/USER-MISC/fix_ttm_mod.cpp b/src/USER-MISC/fix_ttm_mod.cpp index 4bf90fdca..d2fc6ac90 100644 --- a/src/USER-MISC/fix_ttm_mod.cpp +++ b/src/USER-MISC/fix_ttm_mod.cpp @@ -1,922 +1,922 @@ /* ---------------------------------------------------------------------- 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: (in addition to authors of original fix ttm) Sergey Starikov (Joint Institute for High Temperatures of RAS) Vasily Pisarev (Joint Institute for High Temperatures of RAS) ------------------------------------------------------------------------- */ #include "lmptype.h" #include #include #include #include #include "fix_ttm_mod.h" #include "atom.h" #include "force.h" #include "update.h" #include "domain.h" #include "region.h" #include "respa.h" #include "comm.h" #include "random_mars.h" #include "memory.h" #include "error.h" #include "citeme.h" #include "math_const.h" using namespace LAMMPS_NS; using namespace FixConst; using namespace MathConst; #define MAXLINE 1024 static const char cite_fix_ttm_mod[] = "fix ttm/mod command:\n\n" "@article{Pisarev2014,\n" "author = {Pisarev, V. V. and Starikov, S. V.},\n" "title = {{Atomistic simulation of ion track formation in UO2.}},\n" "journal = {J.~Phys.:~Condens.~Matter},\n" "volume = {26},\n" "number = {47},\n" "pages = {475401},\n" "year = {2014}\n" "}\n\n" "@article{Norman2013,\n" "author = {Norman, G. E. and Starikov, S. V. and Stegailov, V. V. and Saitov, I. M. and Zhilyaev, P. A.},\n" "title = {{Atomistic Modeling of Warm Dense Matter in the Two-Temperature State}},\n" "journal = {Contrib.~Plasm.~Phys.},\n" "number = {2},\n" "volume = {53},\n" "pages = {129--139},\n" "year = {2013}\n" "}\n\n"; /* ---------------------------------------------------------------------- */ FixTTMMod::FixTTMMod(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg) { if (lmp->citeme) lmp->citeme->add(cite_fix_ttm_mod); if (narg < 9) error->all(FLERR,"Illegal fix ttm/mod command"); vector_flag = 1; size_vector = 2; global_freq = 1; extvector = 1; nevery = 1; restart_peratom = 1; restart_global = 1; seed = force->inumeric(FLERR,arg[3]); if (seed <= 0) error->all(FLERR,"Invalid random number seed in fix ttm/mod command"); FILE *fpr_2 = force->open_potential(arg[4]); if (fpr_2 == NULL) { char str[128]; sprintf(str,"Cannot open file %s",arg[4]); error->all(FLERR,str); } nxnodes = force->inumeric(FLERR,arg[5]); nynodes = force->inumeric(FLERR,arg[6]); nznodes = force->inumeric(FLERR,arg[7]); if (nxnodes <= 0 || nynodes <= 0 || nznodes <= 0) error->all(FLERR,"Fix ttm/mod number of nodes must be > 0"); FILE *fpr = force->open_potential(arg[8]); if (fpr == NULL) { char str[128]; sprintf(str,"Cannot open file %s",arg[8]); error->all(FLERR,str); } nfileevery = force->inumeric(FLERR,arg[9]); if (nfileevery > 0) { if (narg != 11) error->all(FLERR,"Illegal fix ttm/mod command"); MPI_Comm_rank(world,&me); if (me == 0) { fp = fopen(arg[10],"w"); if (fp == NULL) { char str[128]; sprintf(str,"Cannot open fix ttm/mod file %s",arg[10]); error->one(FLERR,str); } } } char linee[MAXLINE]; double tresh_d; int tresh_i; // C0 (metal) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); esheat_0 = tresh_d; // C1 (metal*10^3) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); esheat_1 = tresh_d; // C2 (metal*10^6) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); esheat_2 = tresh_d; // C3 (metal*10^9) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); esheat_3 = tresh_d; // C4 (metal*10^12) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); esheat_4 = tresh_d; // C_limit fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); C_limit = tresh_d; //Temperature damping factor fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); T_damp = tresh_d; // rho_e fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); electronic_density = tresh_d; //thermal_diffusion fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); el_th_diff = tresh_d; // gamma_p fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); gamma_p = tresh_d; // gamma_s fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); gamma_s = tresh_d; // v0 fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); v_0 = tresh_d; // average intensity of pulse (source of energy) (metal units) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); intensity = tresh_d; // coordinate of 1st surface in x-direction (in box units) - constant fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%d",&tresh_i); surface_l = tresh_i; // coordinate of 2nd surface in x-direction (in box units) - constant fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%d",&tresh_i); surface_r = tresh_i; // skin_layer = intensity is reduced (I=I0*exp[-x/skin_layer]) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%d",&tresh_i); skin_layer = tresh_i; // width of pulse (picoseconds) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); width = tresh_d; // factor of electronic pressure (PF) Pe = PF*Ce*Te fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); pres_factor = tresh_d; // effective free path of electrons (angstrom) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); free_path = tresh_d; // ionic density (ions*angstrom^{-3}) fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); ionic_density = tresh_d; // if movsur = 0: surface is freezed fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%d",&tresh_i); movsur = tresh_i; // electron_temperature_min fgets(linee,MAXLINE,fpr_2); fgets(linee,MAXLINE,fpr_2); sscanf(linee,"%lg",&tresh_d); electron_temperature_min = tresh_d; fclose(fpr_2); //t_surface is determined by electronic temperature (not constant) t_surface_l = surface_l; mult_factor = intensity; duration = 0.0; v_0_sq = v_0*v_0; surface_double = double(t_surface_l)*(domain->xprd/nxnodes); if ((C_limit+esheat_0) < 0.0) error->all(FLERR,"Fix ttm/mod electronic_specific_heat must be >= 0.0"); if (electronic_density <= 0.0) error->all(FLERR,"Fix ttm/mod electronic_density must be > 0.0"); if (gamma_p < 0.0) error->all(FLERR,"Fix ttm/mod gamma_p must be >= 0.0"); if (gamma_s < 0.0) error->all(FLERR,"Fix ttm/mod gamma_s must be >= 0.0"); if (v_0 < 0.0) error->all(FLERR,"Fix ttm/mod v_0 must be >= 0.0"); if (ionic_density <= 0.0) error->all(FLERR,"Fix ttm/mod ionic_density must be > 0.0"); if (surface_l < 0) error->all(FLERR,"Surface coordinates must be >= 0"); if (surface_l >= surface_r) error->all(FLERR, "Left surface coordinate must be less than right surface coordinate"); // initialize Marsaglia RNG with processor-unique seed random = new RanMars(lmp,seed + comm->me); // allocate per-type arrays for force prefactors gfactor1 = new double[atom->ntypes+1]; gfactor2 = new double[atom->ntypes+1]; // allocate 3d grid variables total_nnodes = nxnodes*nynodes*nznodes; memory->create(nsum,nxnodes,nynodes,nznodes,"ttm/mod:nsum"); memory->create(nsum_all,nxnodes,nynodes,nznodes,"ttm/mod:nsum_all"); memory->create(T_initial_set,nxnodes,nynodes,nznodes,"ttm/mod:T_initial_set"); memory->create(sum_vsq,nxnodes,nynodes,nznodes,"ttm/mod:sum_vsq"); memory->create(sum_mass_vsq,nxnodes,nynodes,nznodes,"ttm/mod:sum_mass_vsq"); memory->create(sum_vsq_all,nxnodes,nynodes,nznodes,"ttm/mod:sum_vsq_all"); memory->create(sum_mass_vsq_all,nxnodes,nynodes,nznodes, "ttm/mod:sum_mass_vsq_all"); memory->create(T_electron_old,nxnodes,nynodes,nznodes,"ttm/mod:T_electron_old"); memory->create(T_electron_first,nxnodes,nynodes,nznodes,"ttm/mod:T_electron_first"); memory->create(T_electron,nxnodes,nynodes,nznodes,"ttm/mod:T_electron"); memory->create(net_energy_transfer,nxnodes,nynodes,nznodes, "ttm/mod:net_energy_transfer"); memory->create(net_energy_transfer_all,nxnodes,nynodes,nznodes, "ttm/mod:net_energy_transfer_all"); flangevin = NULL; grow_arrays(atom->nmax); // zero out the flangevin array for (int i = 0; i < atom->nmax; i++) { flangevin[i][0] = 0; flangevin[i][1] = 0; flangevin[i][2] = 0; } atom->add_callback(0); atom->add_callback(1); // set initial electron temperatures from user input file if (me == 0) read_initial_electron_temperatures(fpr); MPI_Bcast(&T_electron[0][0][0],total_nnodes,MPI_DOUBLE,0,world); fclose(fpr); } /* ---------------------------------------------------------------------- */ FixTTMMod::~FixTTMMod() { if (nfileevery && me == 0) fclose(fp); delete random; delete [] gfactor1; delete [] gfactor2; memory->destroy(nsum); memory->destroy(nsum_all); memory->destroy(T_initial_set); memory->destroy(sum_vsq); memory->destroy(sum_mass_vsq); memory->destroy(sum_vsq_all); memory->destroy(sum_mass_vsq_all); memory->destroy(T_electron_first); memory->destroy(T_electron_old); memory->destroy(T_electron); memory->destroy(flangevin); memory->destroy(net_energy_transfer); memory->destroy(net_energy_transfer_all); } /* ---------------------------------------------------------------------- */ int FixTTMMod::setmask() { int mask = 0; mask |= POST_FORCE; mask |= POST_FORCE_RESPA; mask |= END_OF_STEP; return mask; } /* ---------------------------------------------------------------------- */ void FixTTMMod::init() { if (domain->dimension == 2) error->all(FLERR,"Cannot use fix ttm/mod with 2d simulation"); if (domain->nonperiodic != 0) error->all(FLERR,"Cannot use nonperiodic boundares with fix ttm/mod"); if (domain->triclinic) error->all(FLERR,"Cannot use fix ttm/mod with triclinic box"); // set force prefactors for (int i = 1; i <= atom->ntypes; i++) { gfactor1[i] = - gamma_p / force->ftm2v; gfactor2[i] = sqrt(24.0*force->boltz*gamma_p/update->dt/force->mvv2e) / force->ftm2v; } for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) net_energy_transfer_all[ixnode][iynode][iznode] = 0; if (strstr(update->integrate_style,"respa")) nlevels_respa = ((Respa *) update->integrate)->nlevels; } /* ---------------------------------------------------------------------- */ void FixTTMMod::setup(int vflag) { if (strstr(update->integrate_style,"verlet")) post_force_setup(vflag); else { ((Respa *) update->integrate)->copy_flevel_f(nlevels_respa-1); post_force_respa_setup(vflag,nlevels_respa-1,0); ((Respa *) update->integrate)->copy_f_flevel(nlevels_respa-1); } } /* ---------------------------------------------------------------------- */ void FixTTMMod::post_force(int vflag) { double **x = atom->x; double **v = atom->v; double **f = atom->f; int *type = atom->type; int *mask = atom->mask; int nlocal = atom->nlocal; double dx = domain->xprd/nxnodes; double dy = domain->yprd/nynodes; double dz = domain->zprd/nynodes; double gamma1,gamma2; // apply damping and thermostat to all atoms in fix group for (int i = 0; i < nlocal; i++) { if (mask[i] & groupbit) { double xscale = (x[i][0] - domain->boxlo[0])/domain->xprd; double yscale = (x[i][1] - domain->boxlo[1])/domain->yprd; double zscale = (x[i][2] - domain->boxlo[2])/domain->zprd; int ixnode = static_cast(xscale*nxnodes); int iynode = static_cast(yscale*nynodes); int iznode = static_cast(zscale*nznodes); while (ixnode > nxnodes-1) ixnode -= nxnodes; while (iynode > nynodes-1) iynode -= nynodes; while (iznode > nznodes-1) iznode -= nznodes; while (ixnode < 0) ixnode += nxnodes; while (iynode < 0) iynode += nynodes; while (iznode < 0) iznode += nznodes; if (T_electron[ixnode][iynode][iznode] < 0) error->all(FLERR,"Electronic temperature dropped below zero"); double tsqrt = sqrt(T_electron[ixnode][iynode][iznode]); gamma1 = gfactor1[type[i]]; double vsq = v[i][0]*v[i][0] + v[i][1]*v[i][1] + v[i][2]*v[i][2]; if (vsq > v_0_sq) gamma1 *= (gamma_p + gamma_s)/gamma_p; gamma2 = gfactor2[type[i]] * tsqrt; if (ixnode >= surface_l){ if (ixnode < surface_r){ flangevin[i][0] = gamma1*v[i][0] + gamma2*(random->uniform()-0.5); flangevin[i][1] = gamma1*v[i][1] + gamma2*(random->uniform()-0.5); flangevin[i][2] = gamma1*v[i][2] + gamma2*(random->uniform()-0.5); double x_surf = dx*double(surface_l)+dx; double x_at = x[i][0] - domain->boxlo[0]; int right_xnode = ixnode + 1; int right_ynode = iynode + 1; int right_znode = iznode + 1; if (right_xnode == nxnodes) right_xnode = 0; if (right_ynode == nynodes) right_ynode = 0; if (right_znode == nznodes) right_znode = 0; int left_xnode = ixnode - 1; int left_ynode = iynode - 1; int left_znode = iznode - 1; if (left_xnode == -1) left_xnode = nxnodes - 1; if (left_ynode == -1) left_ynode = nynodes - 1; if (left_znode == -1) left_znode = nznodes - 1; double T_i = T_electron[ixnode][iynode][iznode]; double T_ir = T_electron[right_xnode][iynode][iznode]; double T_iu = T_electron[ixnode][right_ynode][iznode]; double T_if = T_electron[ixnode][iynode][right_znode]; double C_i = el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity; double C_ir = el_properties(T_electron[right_xnode][iynode][iznode]).el_heat_capacity; double C_iu = el_properties(T_electron[ixnode][right_ynode][iznode]).el_heat_capacity; double C_if = el_properties(T_electron[ixnode][iynode][right_znode]).el_heat_capacity; double diff_x = (x_at - x_surf)*(x_at - x_surf); diff_x = pow(diff_x,0.5); double len_factor = diff_x/(diff_x+free_path); if (movsur == 1){ if (x_at >= x_surf){ flangevin[i][0] -= pres_factor/ionic_density*((C_ir*T_ir*free_path/(diff_x+free_path)/(diff_x+free_path)) + (len_factor/dx)*(C_ir*T_ir-C_i*T_i)); flangevin[i][1] -= pres_factor/ionic_density/dy*(C_iu*T_iu-C_i*T_i); flangevin[i][2] -= pres_factor/ionic_density/dz*(C_if*T_if-C_i*T_i); } } else{ flangevin[i][0] -= pres_factor/ionic_density/dx*(C_ir*T_ir-C_i*T_i); flangevin[i][1] -= pres_factor/ionic_density/dy*(C_iu*T_iu-C_i*T_i); flangevin[i][2] -= pres_factor/ionic_density/dz*(C_if*T_if-C_i*T_i); } f[i][0] += flangevin[i][0]; f[i][1] += flangevin[i][1]; f[i][2] += flangevin[i][2]; } } if (movsur == 1){ if (ixnode < surface_l){ t_surface_l = ixnode; } } } } MPI_Allreduce(&t_surface_l,&surface_l,1,MPI_INT,MPI_MIN,world); } /* ---------------------------------------------------------------------- */ void FixTTMMod::post_force_setup(int vflag) { double **f = atom->f; int *mask = atom->mask; int nlocal = atom->nlocal; // apply langevin forces that have been stored from previous run for (int i = 0; i < nlocal; i++) { if (mask[i] & groupbit) { f[i][0] += flangevin[i][0]; f[i][1] += flangevin[i][1]; f[i][2] += flangevin[i][2]; } } } /* ---------------------------------------------------------------------- */ void FixTTMMod::post_force_respa(int vflag, int ilevel, int iloop) { if (ilevel == nlevels_respa-1) post_force(vflag); } /* ---------------------------------------------------------------------- */ void FixTTMMod::post_force_respa_setup(int vflag, int ilevel, int iloop) { if (ilevel == nlevels_respa-1) post_force_setup(vflag); } /* ---------------------------------------------------------------------- */ void FixTTMMod::reset_dt() { for (int i = 1; i <= atom->ntypes; i++) gfactor2[i] = sqrt(24.0*force->boltz*gamma_p/update->dt/force->mvv2e) / force->ftm2v; } /* ---------------------------------------------------------------------- read in initial electron temperatures from a user-specified file only called by proc 0 ------------------------------------------------------------------------- */ void FixTTMMod::read_initial_electron_temperatures(FILE *fpr) { char line[MAXLINE]; for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) T_initial_set[ixnode][iynode][iznode] = 0; // read initial electron temperature values from file int ixnode,iynode,iznode; double T_tmp; while (1) { if (fgets(line,MAXLINE,fpr) == NULL) break; sscanf(line,"%d %d %d %lg",&ixnode,&iynode,&iznode,&T_tmp); if (T_tmp < 0.0) error->one(FLERR,"Fix ttm/mod electron temperatures must be >= 0.0"); T_electron[ixnode][iynode][iznode] = T_tmp; T_initial_set[ixnode][iynode][iznode] = 1; } for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) if (T_initial_set[ixnode][iynode][iznode] == 0) error->one(FLERR,"Initial temperatures not all set in fix ttm/mod"); } /* ---------------------------------------------------------------------- */ el_heat_capacity_thermal_conductivity FixTTMMod::el_properties(double T_e) { el_heat_capacity_thermal_conductivity properties; double T_temp = T_e/1000.0, T_reduced = T_damp*T_temp; double T2 = T_temp*T_temp; double T3 = T2*T_temp; double T4 = T3*T_temp; double poly = esheat_0 + esheat_1*T_temp + esheat_2*T2 + esheat_3*T3 + esheat_4*T4; properties.el_heat_capacity = electronic_density*(poly*exp(-T_reduced*T_reduced) + C_limit); // heat capacity properties.el_thermal_conductivity = el_th_diff*properties.el_heat_capacity; // thermal conductivity return properties; } double FixTTMMod::el_sp_heat_integral(double T_e) { double T_temp = T_e/1000.0, T_reduced = T_damp*T_temp; if (T_damp != 0) return electronic_density*(MY_PIS*(3*esheat_4/pow(T_damp,5)+2*esheat_2/pow(T_damp,3)+4*esheat_0/T_damp)*erf(T_reduced)+ 4*esheat_3/pow(T_damp,4)+4*esheat_1/T_damp/T_damp- ((6*esheat_4*T_temp+4*esheat_3)/pow(T_damp,4)+ (4*esheat_1+4*esheat_4*pow(T_temp,3)+4*esheat_3*T_temp*T_temp+4*esheat_2*T_temp)/T_damp/T_damp)*exp(-T_reduced*T_reduced))*125.0+electronic_density*C_limit*T_e; else return electronic_density*((esheat_0 + C_limit)*T_e + esheat_1*T_temp*T_e/2.0 + esheat_2*T_temp*T_temp*T_e/3.0 + esheat_3*pow(T_temp,3)*T_e/4.0 + esheat_4*pow(T_temp,4)*T_e/5.0); } void FixTTMMod::end_of_step() { double **x = atom->x; double **v = atom->v; double *mass = atom->mass; double *rmass = atom->rmass; int *type = atom->type; int *mask = atom->mask; int nlocal = atom->nlocal; if (movsur == 1){ for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++){ double TTT = T_electron[ixnode][iynode][iznode]; if (TTT > 0){ if (ixnode < t_surface_l) t_surface_l = ixnode; } } } for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) net_energy_transfer[ixnode][iynode][iznode] = 0; for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { double xscale = (x[i][0] - domain->boxlo[0])/domain->xprd; double yscale = (x[i][1] - domain->boxlo[1])/domain->yprd; double zscale = (x[i][2] - domain->boxlo[2])/domain->zprd; int ixnode = static_cast(xscale*nxnodes); int iynode = static_cast(yscale*nynodes); int iznode = static_cast(zscale*nznodes); while (ixnode > nxnodes-1) ixnode -= nxnodes; while (iynode > nynodes-1) iynode -= nynodes; while (iznode > nznodes-1) iznode -= nznodes; while (ixnode < 0) ixnode += nxnodes; while (iynode < 0) iynode += nynodes; while (iznode < 0) iznode += nznodes; if (ixnode >= t_surface_l){ if (ixnode < surface_r) net_energy_transfer[ixnode][iynode][iznode] += (flangevin[i][0]*v[i][0] + flangevin[i][1]*v[i][1] + flangevin[i][2]*v[i][2]); } } MPI_Allreduce(&net_energy_transfer[0][0][0], &net_energy_transfer_all[0][0][0], total_nnodes,MPI_DOUBLE,MPI_SUM,world); double dx = domain->xprd/nxnodes; double dy = domain->yprd/nynodes; double dz = domain->zprd/nznodes; double del_vol = dx*dy*dz; double el_specific_heat = 0.0; double el_thermal_conductivity = el_properties(electron_temperature_min).el_thermal_conductivity; for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) { if (el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity > el_thermal_conductivity) el_thermal_conductivity = el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity; if (el_specific_heat > 0.0) { if ((T_electron[ixnode][iynode][iznode] > 0.0) && (el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity < el_specific_heat)) el_specific_heat = el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity; } else if (T_electron[ixnode][iynode][iznode] > 0.0) el_specific_heat = el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity; } // num_inner_timesteps = # of inner steps (thermal solves) // required this MD step to maintain a stable explicit solve int num_inner_timesteps = 1; double inner_dt = update->dt; double stability_criterion = 0.0; for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) T_electron_first[ixnode][iynode][iznode] = T_electron[ixnode][iynode][iznode]; do { for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) T_electron[ixnode][iynode][iznode] = T_electron_first[ixnode][iynode][iznode]; - stability_criterion = 1.0 - - 2.0*inner_dt/el_specific_heat * + stability_criterion = 1.0 - + 2.0*inner_dt/el_specific_heat * + (el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz)); + if (stability_criterion < 0.0) { + inner_dt = 0.25*el_specific_heat / (el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz)); - if (stability_criterion < 0.0) { - inner_dt = 0.25*el_specific_heat / - (el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz)); - } - num_inner_timesteps = static_cast(update->dt/inner_dt) + 1; - inner_dt = update->dt/double(num_inner_timesteps); - if (num_inner_timesteps > 1000000) - error->warning(FLERR,"Too many inner timesteps in fix ttm/mod",0); - for (int ith_inner_timestep = 0; ith_inner_timestep < num_inner_timesteps; - ith_inner_timestep++) { - for (int ixnode = 0; ixnode < nxnodes; ixnode++) - for (int iynode = 0; iynode < nynodes; iynode++) - for (int iznode = 0; iznode < nznodes; iznode++) - T_electron_old[ixnode][iynode][iznode] = - T_electron[ixnode][iynode][iznode]; - // compute new electron T profile - duration = duration + inner_dt; - for (int ixnode = 0; ixnode < nxnodes; ixnode++) - for (int iynode = 0; iynode < nynodes; iynode++) - for (int iznode = 0; iznode < nznodes; iznode++) { - int right_xnode = ixnode + 1; - int right_ynode = iynode + 1; - int right_znode = iznode + 1; - if (right_xnode == nxnodes) right_xnode = 0; - if (right_ynode == nynodes) right_ynode = 0; - if (right_znode == nznodes) right_znode = 0; - int left_xnode = ixnode - 1; - int left_ynode = iynode - 1; - int left_znode = iznode - 1; - if (left_xnode == -1) left_xnode = nxnodes - 1; - if (left_ynode == -1) left_ynode = nynodes - 1; - if (left_znode == -1) left_znode = nznodes - 1; - double skin_layer_d = double(skin_layer); - double ixnode_d = double(ixnode); - double surface_d = double(t_surface_l); - mult_factor = 0.0; - if (duration < width){ - if (ixnode >= t_surface_l) mult_factor = (intensity/(dx*skin_layer_d))*exp((-1.0)*(ixnode_d - surface_d)/skin_layer_d); - } - if (ixnode < t_surface_l) net_energy_transfer_all[ixnode][iynode][iznode] = 0.0; - double cr_vac = 1; - if (T_electron_old[ixnode][iynode][iznode] == 0) cr_vac = 0; - double cr_v_l_x = 1; - if (T_electron_old[left_xnode][iynode][iznode] == 0) cr_v_l_x = 0; - double cr_v_r_x = 1; - if (T_electron_old[right_xnode][iynode][iznode] == 0) cr_v_r_x = 0; - double cr_v_l_y = 1; - if (T_electron_old[ixnode][left_ynode][iznode] == 0) cr_v_l_y = 0; - double cr_v_r_y = 1; - if (T_electron_old[ixnode][right_ynode][iznode] == 0) cr_v_r_y = 0; - double cr_v_l_z = 1; - if (T_electron_old[ixnode][iynode][left_znode] == 0) cr_v_l_z = 0; - double cr_v_r_z = 1; - if (T_electron_old[ixnode][iynode][right_znode] == 0) cr_v_r_z = 0; - if (cr_vac != 0) { - T_electron[ixnode][iynode][iznode] = - T_electron_old[ixnode][iynode][iznode] + - inner_dt/el_properties(T_electron_old[ixnode][iynode][iznode]).el_heat_capacity * - ((cr_v_r_x*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[right_xnode][iynode][iznode]/2.0).el_thermal_conductivity* - (T_electron_old[right_xnode][iynode][iznode]-T_electron_old[ixnode][iynode][iznode])/dx - - cr_v_l_x*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[left_xnode][iynode][iznode]/2.0).el_thermal_conductivity* - (T_electron_old[ixnode][iynode][iznode]-T_electron_old[left_xnode][iynode][iznode])/dx)/dx + - (cr_v_r_y*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][right_ynode][iznode]/2.0).el_thermal_conductivity* - (T_electron_old[ixnode][right_ynode][iznode]-T_electron_old[ixnode][iynode][iznode])/dy - - cr_v_l_y*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][left_ynode][iznode]/2.0).el_thermal_conductivity* - (T_electron_old[ixnode][iynode][iznode]-T_electron_old[ixnode][left_ynode][iznode])/dy)/dy + - (cr_v_r_z*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][iynode][right_znode]/2.0).el_thermal_conductivity* - (T_electron_old[ixnode][iynode][right_znode]-T_electron_old[ixnode][iynode][iznode])/dz - - cr_v_l_z*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][iynode][left_znode]/2.0).el_thermal_conductivity* - (T_electron_old[ixnode][iynode][iznode]-T_electron_old[ixnode][iynode][left_znode])/dz)/dz); - T_electron[ixnode][iynode][iznode]+=inner_dt/el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity* - (mult_factor - - net_energy_transfer_all[ixnode][iynode][iznode]/del_vol); - } - else T_electron[ixnode][iynode][iznode] = - T_electron_old[ixnode][iynode][iznode]; - if ((T_electron[ixnode][iynode][iznode] > 0.0) && (T_electron[ixnode][iynode][iznode] < electron_temperature_min)) - T_electron[ixnode][iynode][iznode] = T_electron[ixnode][iynode][iznode] + 0.5*(electron_temperature_min - T_electron[ixnode][iynode][iznode]); - - if (el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity > el_thermal_conductivity) - el_thermal_conductivity = el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity; - if ((T_electron[ixnode][iynode][iznode] > 0.0) && (el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity < el_specific_heat)) - el_specific_heat = el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity; + } + num_inner_timesteps = static_cast(update->dt/inner_dt) + 1; + inner_dt = update->dt/double(num_inner_timesteps); + if (num_inner_timesteps > 1000000) + error->warning(FLERR,"Too many inner timesteps in fix ttm/mod",0); + for (int ith_inner_timestep = 0; ith_inner_timestep < num_inner_timesteps; + ith_inner_timestep++) { + for (int ixnode = 0; ixnode < nxnodes; ixnode++) + for (int iynode = 0; iynode < nynodes; iynode++) + for (int iznode = 0; iznode < nznodes; iznode++) + T_electron_old[ixnode][iynode][iznode] = + T_electron[ixnode][iynode][iznode]; + // compute new electron T profile + duration = duration + inner_dt; + for (int ixnode = 0; ixnode < nxnodes; ixnode++) + for (int iynode = 0; iynode < nynodes; iynode++) + for (int iznode = 0; iznode < nznodes; iznode++) { + int right_xnode = ixnode + 1; + int right_ynode = iynode + 1; + int right_znode = iznode + 1; + if (right_xnode == nxnodes) right_xnode = 0; + if (right_ynode == nynodes) right_ynode = 0; + if (right_znode == nznodes) right_znode = 0; + int left_xnode = ixnode - 1; + int left_ynode = iynode - 1; + int left_znode = iznode - 1; + if (left_xnode == -1) left_xnode = nxnodes - 1; + if (left_ynode == -1) left_ynode = nynodes - 1; + if (left_znode == -1) left_znode = nznodes - 1; + double skin_layer_d = double(skin_layer); + double ixnode_d = double(ixnode); + double surface_d = double(t_surface_l); + mult_factor = 0.0; + if (duration < width){ + if (ixnode >= t_surface_l) mult_factor = (intensity/(dx*skin_layer_d))*exp((-1.0)*(ixnode_d - surface_d)/skin_layer_d); } - } - stability_criterion = 1.0 - - 2.0*inner_dt/el_specific_heat * - (el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz)); + if (ixnode < t_surface_l) net_energy_transfer_all[ixnode][iynode][iznode] = 0.0; + double cr_vac = 1; + if (T_electron_old[ixnode][iynode][iznode] == 0) cr_vac = 0; + double cr_v_l_x = 1; + if (T_electron_old[left_xnode][iynode][iznode] == 0) cr_v_l_x = 0; + double cr_v_r_x = 1; + if (T_electron_old[right_xnode][iynode][iznode] == 0) cr_v_r_x = 0; + double cr_v_l_y = 1; + if (T_electron_old[ixnode][left_ynode][iznode] == 0) cr_v_l_y = 0; + double cr_v_r_y = 1; + if (T_electron_old[ixnode][right_ynode][iznode] == 0) cr_v_r_y = 0; + double cr_v_l_z = 1; + if (T_electron_old[ixnode][iynode][left_znode] == 0) cr_v_l_z = 0; + double cr_v_r_z = 1; + if (T_electron_old[ixnode][iynode][right_znode] == 0) cr_v_r_z = 0; + if (cr_vac != 0) { + T_electron[ixnode][iynode][iznode] = + T_electron_old[ixnode][iynode][iznode] + + inner_dt/el_properties(T_electron_old[ixnode][iynode][iznode]).el_heat_capacity * + ((cr_v_r_x*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[right_xnode][iynode][iznode]/2.0).el_thermal_conductivity* + (T_electron_old[right_xnode][iynode][iznode]-T_electron_old[ixnode][iynode][iznode])/dx - + cr_v_l_x*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[left_xnode][iynode][iznode]/2.0).el_thermal_conductivity* + (T_electron_old[ixnode][iynode][iznode]-T_electron_old[left_xnode][iynode][iznode])/dx)/dx + + (cr_v_r_y*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][right_ynode][iznode]/2.0).el_thermal_conductivity* + (T_electron_old[ixnode][right_ynode][iznode]-T_electron_old[ixnode][iynode][iznode])/dy - + cr_v_l_y*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][left_ynode][iznode]/2.0).el_thermal_conductivity* + (T_electron_old[ixnode][iynode][iznode]-T_electron_old[ixnode][left_ynode][iznode])/dy)/dy + + (cr_v_r_z*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][iynode][right_znode]/2.0).el_thermal_conductivity* + (T_electron_old[ixnode][iynode][right_znode]-T_electron_old[ixnode][iynode][iznode])/dz - + cr_v_l_z*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][iynode][left_znode]/2.0).el_thermal_conductivity* + (T_electron_old[ixnode][iynode][iznode]-T_electron_old[ixnode][iynode][left_znode])/dz)/dz); + T_electron[ixnode][iynode][iznode]+=inner_dt/el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity* + (mult_factor - + net_energy_transfer_all[ixnode][iynode][iznode]/del_vol); + } + else T_electron[ixnode][iynode][iznode] = + T_electron_old[ixnode][iynode][iznode]; + if ((T_electron[ixnode][iynode][iznode] > 0.0) && (T_electron[ixnode][iynode][iznode] < electron_temperature_min)) + T_electron[ixnode][iynode][iznode] = T_electron[ixnode][iynode][iznode] + 0.5*(electron_temperature_min - T_electron[ixnode][iynode][iznode]); + + if (el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity > el_thermal_conductivity) + el_thermal_conductivity = el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity; + if ((T_electron[ixnode][iynode][iznode] > 0.0) && (el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity < el_specific_heat)) + el_specific_heat = el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity; + } + } + stability_criterion = 1.0 - + 2.0*inner_dt/el_specific_heat * + (el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz)); - } while (stability_criterion < 0.0); + } while (stability_criterion < 0.0); // output nodal temperatures for current timestep if ((nfileevery) && !(update->ntimestep % nfileevery)) { // compute atomic Ta for each grid point for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) { nsum[ixnode][iynode][iznode] = 0; nsum_all[ixnode][iynode][iznode] = 0; sum_vsq[ixnode][iynode][iznode] = 0.0; sum_mass_vsq[ixnode][iynode][iznode] = 0.0; sum_vsq_all[ixnode][iynode][iznode] = 0.0; sum_mass_vsq_all[ixnode][iynode][iznode] = 0.0; } double massone; for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { if (rmass) massone = rmass[i]; else massone = mass[type[i]]; double xscale = (x[i][0] - domain->boxlo[0])/domain->xprd; double yscale = (x[i][1] - domain->boxlo[1])/domain->yprd; double zscale = (x[i][2] - domain->boxlo[2])/domain->zprd; int ixnode = static_cast(xscale*nxnodes); int iynode = static_cast(yscale*nynodes); int iznode = static_cast(zscale*nznodes); while (ixnode > nxnodes-1) ixnode -= nxnodes; while (iynode > nynodes-1) iynode -= nynodes; while (iznode > nznodes-1) iznode -= nznodes; while (ixnode < 0) ixnode += nxnodes; while (iynode < 0) iynode += nynodes; while (iznode < 0) iznode += nznodes; double vsq = v[i][0]*v[i][0] + v[i][1]*v[i][1] + v[i][2]*v[i][2]; nsum[ixnode][iynode][iznode] += 1; sum_vsq[ixnode][iynode][iznode] += vsq; sum_mass_vsq[ixnode][iynode][iznode] += massone*vsq; } MPI_Allreduce(&nsum[0][0][0],&nsum_all[0][0][0],total_nnodes, MPI_INT,MPI_SUM,world); MPI_Allreduce(&sum_vsq[0][0][0],&sum_vsq_all[0][0][0],total_nnodes, MPI_DOUBLE,MPI_SUM,world); MPI_Allreduce(&sum_mass_vsq[0][0][0],&sum_mass_vsq_all[0][0][0], total_nnodes,MPI_DOUBLE,MPI_SUM,world); MPI_Allreduce(&t_surface_l,&surface_l, 1,MPI_INT,MPI_MIN,world); if (me == 0) { fprintf(fp,BIGINT_FORMAT,update->ntimestep); double T_a; for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) { T_a = 0; if (nsum_all[ixnode][iynode][iznode] > 0){ T_a = sum_mass_vsq_all[ixnode][iynode][iznode]/ (3.0*force->boltz*nsum_all[ixnode][iynode][iznode]/force->mvv2e); if (movsur == 1){ if (T_electron[ixnode][iynode][iznode]==0.0) T_electron[ixnode][iynode][iznode] = electron_temperature_min; } } fprintf(fp," %f",T_a); } fprintf(fp,"\t"); for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) fprintf(fp,"%f ",T_electron[ixnode][iynode][iznode]); fprintf(fp,"\n"); } } } /* ---------------------------------------------------------------------- memory usage of 3d grid ------------------------------------------------------------------------- */ double FixTTMMod::memory_usage() { double bytes = 0.0; bytes += 5*total_nnodes * sizeof(int); bytes += 14*total_nnodes * sizeof(double); return bytes; } /* ---------------------------------------------------------------------- */ void FixTTMMod::grow_arrays(int ngrow) { memory->grow(flangevin,ngrow,3,"ttm/mod:flangevin"); } /* ---------------------------------------------------------------------- return the energy of the electronic subsystem or the net_energy transfer between the subsystems ------------------------------------------------------------------------- */ double FixTTMMod::compute_vector(int n) { double e_energy = 0.0; double transfer_energy = 0.0; double dx = domain->xprd/nxnodes; double dy = domain->yprd/nynodes; double dz = domain->zprd/nznodes; double del_vol = dx*dy*dz; for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) { e_energy += el_sp_heat_integral(T_electron[ixnode][iynode][iznode])*del_vol; transfer_energy += net_energy_transfer_all[ixnode][iynode][iznode]*update->dt; } if (n == 0) return e_energy; if (n == 1) return transfer_energy; return 0.0; } /* ---------------------------------------------------------------------- pack entire state of Fix into one write ------------------------------------------------------------------------- */ void FixTTMMod::write_restart(FILE *fp) { double *rlist; memory->create(rlist,nxnodes*nynodes*nznodes+1,"ttm/mod:rlist"); int n = 0; rlist[n++] = seed; for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) rlist[n++] = T_electron[ixnode][iynode][iznode]; if (comm->me == 0) { int size = n * sizeof(double); fwrite(&size,sizeof(int),1,fp); fwrite(rlist,sizeof(double),n,fp); } memory->destroy(rlist); } /* ---------------------------------------------------------------------- use state info from restart file to restart the Fix ------------------------------------------------------------------------- */ void FixTTMMod::restart(char *buf) { int n = 0; double *rlist = (double *) buf; // the seed must be changed from the initial seed seed = static_cast (0.5*rlist[n++]); for (int ixnode = 0; ixnode < nxnodes; ixnode++) for (int iynode = 0; iynode < nynodes; iynode++) for (int iznode = 0; iznode < nznodes; iznode++) T_electron[ixnode][iynode][iznode] = rlist[n++]; delete random; random = new RanMars(lmp,seed+comm->me); } /* ---------------------------------------------------------------------- pack values in local atom-based arrays for restart file ------------------------------------------------------------------------- */ int FixTTMMod::pack_restart(int i, double *buf) { buf[0] = 4; buf[1] = flangevin[i][0]; buf[2] = flangevin[i][1]; buf[3] = flangevin[i][2]; return 4; } /* ---------------------------------------------------------------------- unpack values from atom->extra array to restart the fix ------------------------------------------------------------------------- */ void FixTTMMod::unpack_restart(int nlocal, int nth) { 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++; flangevin[nlocal][0] = extra[nlocal][m++]; flangevin[nlocal][1] = extra[nlocal][m++]; flangevin[nlocal][2] = extra[nlocal][m++]; } /* ---------------------------------------------------------------------- maxsize of any atom's restart data ------------------------------------------------------------------------- */ int FixTTMMod::maxsize_restart() { return 4; } /* ---------------------------------------------------------------------- size of atom nlocal's restart data ------------------------------------------------------------------------- */ int FixTTMMod::size_restart(int nlocal) { return 4; }