Merge pull request #2349 from akohlmey/more-progguide-updates

More Programmer docs updates and related code refactoring
This commit is contained in:
Axel Kohlmeyer 2020-09-09 14:34:41 -04:00 committed by GitHub
commit cdd9d693ad
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
37 changed files with 2638 additions and 1706 deletions

View File

@ -417,7 +417,6 @@ INPUT = @LAMMPS_SOURCE_DIR@/utils.cpp \
@LAMMPS_SOURCE_DIR@/lammps.cpp \
@LAMMPS_SOURCE_DIR@/lammps.h \
@LAMMPS_SOURCE_DIR@/lmptype.h \
@LAMMPS_SOURCE_DIR@/pointers.h \
@LAMMPS_SOURCE_DIR@/atom.cpp \
@LAMMPS_SOURCE_DIR@/atom.h \
@LAMMPS_SOURCE_DIR@/input.cpp \
@ -428,6 +427,10 @@ INPUT = @LAMMPS_SOURCE_DIR@/utils.cpp \
@LAMMPS_SOURCE_DIR@/text_file_reader.h \
@LAMMPS_SOURCE_DIR@/potential_file_reader.cpp \
@LAMMPS_SOURCE_DIR@/potential_file_reader.h \
@LAMMPS_SOURCE_DIR@/my_page.cpp \
@LAMMPS_SOURCE_DIR@/my_page.h \
@LAMMPS_SOURCE_DIR@/my_pool_chunk.cpp \
@LAMMPS_SOURCE_DIR@/my_pool_chunk.h \
# The EXCLUDE_SYMLINKS tag can be used to select whether or not files or
# directories that are symbolic links (a Unix file system feature) are excluded

View File

@ -72,8 +72,6 @@ every LAMMPS command.
pg_library
Modify
pg_developer
.. pg_modify
.. pg_base
.. toctree::
:caption: Index

9
doc/src/pg_atom.rst Normal file
View File

@ -0,0 +1,9 @@
LAMMPS Atom and AtomVec Base Classes
************************************
.. doxygenclass:: LAMMPS_NS::Atom
:project: progguide
:members:

View File

@ -0,0 +1,38 @@
LAMMPS C++ base classes
=======================
LAMMPS is designed to be used as a C++ class library where one can set
up and drive a simulation through creating a class instance and then
calling some abstract operations or commands on that class or its member
class instances. These are interfaced to the :doc:`C library API
<pg_library>`, which providing an additional level of abstraction
simplification for common operations. The C API is also the basis for
calling LAMMPS from Python or Fortran.
When used from a C++ program, most of the symbols and functions in
LAMMPS are wrapped into the ``LAMMPS_NS`` namespace so they will not
collide with your own classes or other libraries. This, however, does
not extend to the additional libraries bundled with LAMMPS in the lib
folder and some of the low-level code of some packages.
Behind the scenes this is implemented through inheritance and
polymorphism where base classes define the abstract interface and
derived classes provide the specialized implementation for specific
models or optimizations or ports to accelerator platforms. This
document will provide an outline of the fundamental class hierarchy and
some selected examples for derived classes of specific models.
.. note::
Please see the :ref:`note about thread-safety <thread-safety>`
in the library Howto doc page.
-----------------------------------
.. toctree::
:caption: Individual Base Classes
:name: lammpsbase
pg_lammps
pg_atom
pg_input

236
doc/src/pg_dev_flow.rst Normal file
View File

@ -0,0 +1,236 @@
How a timestep works
====================
The first and most fundamental operation within LAMMPS to understand is
how a timestep is structured. Timestepping is performed by calling
methods of the Integrate class instance within the Update class. Since
Integrate is a base class, it will point to an instance of a derived
class corresponding to what is selected by the :doc:`run_style
<run_style>` input script command.
In this section, the timestep implemented by the Verlet class is
described. A similar timestep protocol is implemented by the Respa
class, for the r-RESPA hierarchical timestepping method.
The Min base class performs energy minimization, so does not perform a
literal timestep. But it has logic similar to what is described here,
to compute forces and invoke fixes at each iteration of a minimization.
Differences between time integration and minimization are highlighted at
the end of this section.
The Verlet class is encoded in the ``src/verlet.cpp`` and ``verlet.h``
files. It implements the velocity-Verlet timestepping algorithm. The
workhorse method is ``Verlet::run()``, but first we highlight several
other methods in the class.
- The ``init()`` method is called at the beginning of each dynamics
run. It simply sets some internal flags, based on user settings in
other parts of the code.
- The ``setup()`` or ``setup_minimal()`` methods are also called before
each run. The velocity-Verlet method requires current forces be
calculated before the first timestep, so these routines compute
forces due to all atomic interactions, using the same logic that
appears in the timestepping described next. A few fixes are also
invoked, using the mechanism described in the next section. Various
counters are also initialized before the run begins. The
``setup_minimal()`` method is a variant that has a flag for performing
less setup. This is used when runs are continued and information
from the previous run is still valid. For example, if repeated
short LAMMPS runs are being invoked, interleaved by other commands,
via the *pre no* and *every* options of the run command, the
``setup_minimal()`` method is used.
- The ``force_clear()`` method initializes force and other arrays to
zero before each timestep, so that forces (torques, etc) can be
accumulated.
Now for the ``Verlet::run()`` method. Its basic structure in hi-level pseudo
code is shown below. In the actual code in ``src/verlet.cpp`` some of
these operations are conditionally invoked.
.. code-block:: python
loop over N timesteps:
if timeout condition: break
ev_set()
fix->initial_integrate()
fix->post_integrate()
nflag = neighbor->decide()
if nflag:
fix->pre_exchange()
domain->pbc()
domain->reset_box()
comm->setup()
neighbor->setup_bins()
comm->exchange()
comm->borders()
fix->pre_neighbor()
neighbor->build()
fix->post_neighbor()
else:
comm->forward_comm()
force_clear()
fix->pre_force()
pair->compute()
bond->compute()
angle->compute()
dihedral->compute()
improper->compute()
kspace->compute()
fix->pre_reverse()
comm->reverse_comm()
fix->post_force()
fix->final_integrate()
fix->end_of_step()
if any output on this step:
output->write()
# after loop
fix->post_run()
The ``ev_set()`` method (in the parent Integrate class), sets two flags
(*eflag* and *vflag*) for energy and virial computation. Each flag
encodes whether global and/or per-atom energy and virial should be
calculated on this timestep, because some fix or variable or output will
need it. These flags are passed to the various methods that compute
particle interactions, so that they either compute and tally the
corresponding data or can skip the extra calculations if the energy and
virial are not needed. See the comments for the ``Integrate::ev_set()``
method which document the flag values.
At various points of the timestep, fixes are invoked,
e.g. ``fix->initial_integrate()``. In the code, this is actually done
via the Modify class which stores all the Fix objects and lists of which
should be invoked at what point in the timestep. Fixes are the LAMMPS
mechanism for tailoring the operations of a timestep for a particular
simulation. As described elsewhere, each fix has one or more methods,
each of which is invoked at a specific stage of the timestep, as show in
the timestep pseudo-code. All the active fixes defined in an input
script, that are flagged to have an ``initial_integrate()`` method are
invoked at the beginning of each timestep. Examples are :doc:`fix nve
<fix_nve>` or :doc:`fix nvt or fix npt <fix_nh>` which perform the
start-of-timestep velocity-Verlet integration operations to update
velocities by a half-step, and coordinates by a full step. The
``post_integrate()`` method is next for operations that need to happen
immediately after those updates. Only a few fixes use this, e.g. to
reflect particles off box boundaries in the :doc:`FixWallReflect class
<fix_wall_reflect>`.
The ``decide()`` method in the Neighbor class determines whether
neighbor lists need to be rebuilt on the current timestep (conditions
can be changed using the :doc:`neigh_modify every/delay/check
<neigh_modify>` command. If not, coordinates of ghost atoms are
acquired by each processor via the ``forward_comm()`` method of the Comm
class. If neighbor lists need to be built, several operations within
the inner if clause of the pseudo-code are first invoked. The
``pre_exchange()`` method of any defined fixes is invoked first.
Typically this inserts or deletes particles from the system.
Periodic boundary conditions are then applied by the Domain class via
its ``pbc()`` method to remap particles that have moved outside the
simulation box back into the box. Note that this is not done every
timestep, but only when neighbor lists are rebuilt. This is so that
each processor's sub-domain will have consistent (nearby) atom
coordinates for its owned and ghost atoms. It is also why dumped atom
coordinates may be slightly outside the simulation box if not dumped
on a step where the neighbor lists are rebuilt.
The box boundaries are then reset (if needed) via the ``reset_box()``
method of the Domain class, e.g. if box boundaries are shrink-wrapped to
current particle coordinates. A change in the box size or shape
requires internal information for communicating ghost atoms (Comm class)
and neighbor list bins (Neighbor class) be updated. The ``setup()``
method of the Comm class and ``setup_bins()`` method of the Neighbor
class perform the update.
The code is now ready to migrate atoms that have left a processor's
geometric sub-domain to new processors. The ``exchange()`` method of
the Comm class performs this operation. The ``borders()`` method of the
Comm class then identifies ghost atoms surrounding each processor's
sub-domain and communicates ghost atom information to neighboring
processors. It does this by looping over all the atoms owned by a
processor to make lists of those to send to each neighbor processor. On
subsequent timesteps, the lists are used by the ``Comm::forward_comm()``
method.
Fixes with a ``pre_neighbor()`` method are then called. These typically
re-build some data structure stored by the fix that depends on the
current atoms owned by each processor.
Now that each processor has a current list of its owned and ghost
atoms, LAMMPS is ready to rebuild neighbor lists via the ``build()``
method of the Neighbor class. This is typically done by binning all
owned and ghost atoms, and scanning a stencil of bins around each
owned atom's bin to make a Verlet list of neighboring atoms within the
force cutoff plus neighbor skin distance.
In the next portion of the timestep, all interaction forces between
particles are computed, after zeroing the per-atom force vector via the
``force_clear()`` method. If the newton flag is set to *on* by the
newton command, forces are added to both owned and ghost atoms, otherwise
only to owned (aka local) atoms.
Pairwise forces are calculated first, which enables the global virial
(if requested) to be calculated cheaply (at O(N) cost instead of O(N**2)
at the end of the ``Pair::compute()`` method), by a dot product of atom
coordinates and forces. By including owned and ghost atoms in the dot
product, the effect of periodic boundary conditions is correctly
accounted for. Molecular topology interactions (bonds, angles,
dihedrals, impropers) are calculated next (if supported by the current
atom style). The final contribution is from long-range Coulombic
interactions, invoked by the KSpace class.
The ``pre_reverse()`` method in fixes is used for operations that have to
be done *before* the upcoming reverse communication (e.g. to perform
additional data transfers or reductions for data computed during the
force computation and stored with ghost atoms).
If the newton flag is on, forces on ghost atoms are communicated and
summed back to their corresponding owned atoms. The ``reverse_comm()``
method of the Comm class performs this operation, which is essentially
the inverse operation of sending copies of owned atom coordinates to
other processor's ghost atoms.
At this point in the timestep, the total force on each (local) atom is
known. Additional force constraints (external forces, SHAKE, etc) are
applied by Fixes that have a ``post_force()`` method. The second half
of the velocity-Verlet integration, ``final_integrate()`` is then
performed (another half-step update of the velocities) via fixes like
nve, nvt, npt.
At the end of the timestep, fixes that contain an ``end_of_step()``
method are invoked. These typically perform a diagnostic calculation,
e.g. the ave/time and ave/spatial fixes. The final operation of the
timestep is to perform any requested output, via the ``write()`` method
of the Output class. There are 3 kinds of LAMMPS output: thermodynamic
output to the screen and log file, snapshots of atom data to a dump
file, and restart files. See the :doc:`thermo_style <thermo_style>`,
:doc:`dump <dump>`, and :doc:`restart <restart>` commands for more
details.
The the flow of control during energy minimization iterations is
similar to that of a molecular dynamics timestep. Forces are computed,
neighbor lists are built as needed, atoms migrate to new processors, and
atom coordinates and forces are communicated to neighboring processors.
The only difference is what Fix class operations are invoked when. Only
a subset of LAMMPS fixes are useful during energy minimization, as
explained in their individual doc pages. The relevant Fix class methods
are ``min_pre_exchange()``, ``min_pre_force()``, and ``min_post_force()``.
Each fix is invoked at the appropriate place within the minimization
iteration. For example, the ``min_post_force()`` method is analogous to
the ``post_force()`` method for dynamics; it is used to alter or constrain
forces on each atom, which affects the minimization procedure.
After all iterations are completed there is a ``cleanup`` step which
calls the ``post_run()`` method of fixes to perform operations only required
at the end of a calculations (like freeing temporary storage or creating
final outputs).

250
doc/src/pg_dev_org.rst Normal file
View File

@ -0,0 +1,250 @@
LAMMPS source files
===================
The source files of the LAMMPS code are found in two
directories of the distribution: ``src`` and ``lib``.
Most of the code is C++ but there are small numbers of files
in several other languages.
The core of the code is located in the
``src`` folder and its sub-directories.
A sizable number of these files are in the ``src`` directory
itself, but there are plenty of :doc:`packages <Packages>`, which can be
included or excluded when LAMMPS is built. See the :doc:`Include
packages in build <Build_package>` section of the manual for more
information about that part of the build process. LAMMPS currently
supports building with :doc:`conventional makefiles <Build_make>` and
through :doc:`CMake <Build_cmake>` which differ in how packages are
enabled or disabled for a LAMMPS binary. The source files for each
package are in all-uppercase sub-directories of the ``src`` folder, for
example ``src/MOLECULE`` or ``src/USER-MISC``. The ``src/STUBS``
sub-directory is not a package but contains a dummy MPI library, that is
used when building a serial version of the code. The ``src/MAKE``
directory contains makefiles with settings and flags for a variety of
configuration and machines for the build process with traditional
makefiles.
The ``lib`` directory contains the source code for several supporting
libraries or files with configuration settings to use globally installed
libraries, that are required by some of the optional packages.
Each sub-directory, like ``lib/poems`` or ``lib/gpu``, contains the
source files, some of which are in different languages such as Fortran
or CUDA. These libraries are linked to during a LAMMPS build, if the
corresponding package is installed.
LAMMPS C++ source files almost always come in pairs, such as
``src/run.cpp`` (implementation file) and ``src/run.h`` (header file).
Each pair of files defines a C++
class, for example the :cpp:class:`LAMMPS_NS::Run` class which contains
the code invoked by the :doc:`run <run>` command in a LAMMPS input script.
As this example illustrates, source file and class names often have a
one-to-one correspondence with a command used in a LAMMPS input script.
Some source files and classes do not have a corresponding input script
command, e.g. ``src/force.cpp`` and the :cpp:class:`LAMMPS_NS::Force`
class. They are discussed in the next section.
A small number of C++ classes and utility functions are implemented with
only a ``.h`` file. Examples are the Pointer class or the MathVec functions.
LAMMPS class topology
=====================
Though LAMMPS has a lot of source files and classes, its class topology
is relative flat, as outlined in the :ref:`class-topology` figure. Each
name refers to a class and has a pair of associated source files in the
``src`` folder, for example the class :cpp:class:`LAMMPS_NS::Memory`
corresponds to the files ``memory.cpp`` and ``memory.h``, or the class
:cpp:class:`LAMMPS_NS::AtomVec` corresponds to the files
``atom_vec.cpp`` and ``atom_vec.h``. Full lines in the figure represent
compositing: that is the class to the left holds a pointer to an
instance of the class to the right. Dashed lines instead represent
inheritance: the class to the right is derived from the class on the
left. Classes with a red boundary are not instantiated directly, but
they represent the base classes for "styles". Those "styles" make up
the bulk of the LAMMPS code and only a few typical examples are included
in the figure for demonstration purposes.
.. _class-topology:
.. figure:: JPG/lammps-classes.png
LAMMPS class topology
This figure shows some of the relations of the base classes of the
LAMMPS simulation package. Full lines indicate that a class holds an
instance of the class it is pointing to; dashed lines point to
derived classes that are given as examples of what classes may be
instantiated during a LAMMPS run based on the input commands and
accessed through the API define by their respective base classes. At
the core is the :cpp:class:`LAMMPS <LAMMPS_NS::LAMMPS>` class, which
holds pointers to class instances with specific purposes. Those may
hold instances of other classes, sometimes directly, or only
temporarily, sometimes as derived classes or derived classes or
derived classes, which may also hold instances of other classes.
The :cpp:class:`LAMMPS_NS::LAMMPS` class is the topmost class and
represents what is referred to an "instance" of LAMMPS. It is a
composite holding references to instances of other core classes
providing the core functionality of the MD engine in LAMMPS and through
them abstractions of the required operations. The constructor of the
LAMMPS class will instantiate those instances, process the command line
flags, initialize MPI (if not already done) and set up file pointers for
input and output. The destructor will shut everything down and free all
associated memory. Thus code for the standalone LAMMPS executable in
``main.cpp`` simply initializes MPI, instantiates a single instance of
LAMMPS, and passes it the command line flags and input script. It
deletes the LAMMPS instance after the method reading the input returns
and shuts down the MPI environment before it exits the executable.
The :cpp:class:`LAMMPS_NS::Pointers` is not shown in the
:ref:`class-topology` figure, it holds references to members of the
`LAMMPS_NS::LAMMPS`, so that all classes derived from
:cpp:class:`LAMMPS_NS::Pointers` have direct access to those reference.
From the class topology all classes with blue boundary are referenced in
this class and all classes in the second and third columns, that are not
listed as derived classes are instead derived from
:cpp:class:`LAMMPS_NS::Pointers`.
Since all storage is encapsulated, the LAMMPS class can also be
instantiated multiple times by a calling code, and that can be either
simultaneously or consecutively. When running in parallel with MPI,
care has to be taken, that suitable communicators are used to not
create conflicts between different instances.
The LAMMPS class currently holds instances of 19 classes representing
different core functionalities There are a handful of virtual parent
classes in LAMMPS that define what LAMMPS calls ``styles``. They are
shaded red in the :ref:`class-topology` figure. Each of these are
parents of a number of child classes that implement the interface
defined by the parent class. There are two main categories of these
``styles``: some may only have one instance active at a time (e.g. atom,
pair, bond, angle, dihedral, improper, kspace, comm) and there is a
dedicated pointer variable in the composite class that manages them.
Setups that require a mix of different such styles have to use a
*hybrid* class that manages and forwards calls to the corresponding
sub-styles for the designated subset of atoms or data. or the composite
class may have lists of class instances, e.g. Modify handles lists of
compute and fix styles, while Output handles dumps class instances.
The exception to this scheme are the ``command`` style classes. These
implement specific commands that can be invoked before, after, or between
runs or are commands which launch a simulation. For these an instance
of the class is created, its command() method called and then, after
completion, the class instance deleted. Examples for this are the
create_box, create_atoms, minimize, run, or velocity command styles.
For all those ``styles`` certain naming conventions are employed: for
the fix nve command the class is called FixNVE and the files are
``fix_nve.h`` and ``fix_nve.cpp``. Similarly for fix ave/time we have
FixAveTime and ``fix_ave_time.h`` and ``fix_ave_time.cpp``. Style names
are lower case and without spaces or special characters. A suffix or
multiple appended with a forward slash '/' denotes a variant of the
corresponding class without the suffix. To connect the style name and
the class name, LAMMPS uses macros like the following ATOM\_CLASS,
PAIR\_CLASS, BOND\_CLASS, REGION\_CLASS, FIX\_CLASS, COMPUTE\_CLASS,
or DUMP\_CLASS in the corresponding header file. During compilation
files with the pattern ``style_name.h`` are created that contain include
statements including all headers of all styles of a given type that
are currently active (or "installed).
More details on individual classes in the :ref:`class-topology` are as
follows:
- The Memory class handles allocation of all large vectors and arrays.
- The Error class prints all error and warning messages.
- The Universe class sets up partitions of processors so that multiple
simulations can be run, each on a subset of the processors allocated
for a run, e.g. by the mpirun command.
- The Input class reads and processes input input strings and files,
stores variables, and invokes :doc:`commands <Commands_all>`.
- As discussed above, command style classes are directly derived from
the Pointers class. They provide input script commands that perform
one-time operations before/after/between simulations or which invoke a
simulation. They are instantiated from within the Input class,
invoked, then immediately destructed.
- The Finish class is instantiated to print statistics to the screen
after a simulation is performed, by commands like run and minimize.
- The Special class walks the bond topology of a molecular system to
find first, second, third neighbors of each atom. It is invoked by
several commands, like :doc:`read_data <read_data>`,
:doc:`read_restart <read_restart>`, or :doc:`replicate <replicate>`.
- The Atom class stores per-atom properties associated with atom styles.
More precisely, they are allocated and managed by a class derived from
the AtomVec class, and the Atom class simply stores pointers to them.
The classes derived from AtomVec represent the different atom styles
and they are instantiated through the :doc:`atom_style <atom_style>`
command.
- The Update class holds instances of an integrator and a minimizer
class. The Integrate class is a parent style for the Verlet and
r-RESPA time integrators, as defined by the :doc:`run_style
<run_style>` command. The Min class is a parent style for various
energy minimizers.
- The Neighbor class builds and stores neighbor lists. The NeighList
class stores a single list (for all atoms). A NeighRequest class
instance is created by pair, fix, or compute styles when they need a
particular kind of neighbor list and use the NeighRequest properties
to select the neighbor list settings for the given request. There can
be multiple instances of the NeighRequest class and the Neighbor class
will try to optimize how they are computed by creating copies or
sub-lists where possible.
- The Comm class performs inter-processor communication, typically of
ghost atom information. This usually involves MPI message exchanges
with 6 neighboring processors in the 3d logical grid of processors
mapped to the simulation box. There are two :doc:`communication styles
<comm_style>` enabling different ways to do the domain decomposition.
Sometimes the Irregular class is used, when atoms may migrate to
arbitrary processors.
- The Domain class stores the simulation box geometry, as well as
geometric Regions and any user definition of a Lattice. The latter
are defined by the :doc:`region <region>` and :doc:`lattice <lattice>`
commands in an input script.
- The Force class computes various forces between atoms. The Pair
parent class is for non-bonded or pair-wise forces, which in LAMMPS
also includes many-body forces such as the Tersoff 3-body potential if
those are computed by walking pairwise neighbor lists. The Bond,
Angle, Dihedral, Improper parent classes are styles for bonded
interactions within a static molecular topology. The KSpace parent
class is for computing long-range Coulombic interactions. One of its
child classes, PPPM, uses the FFT3D and Remap classes to redistribute
and communicate grid-based information across the parallel processors.
- The Modify class stores lists of class instances derived from the
:doc:`Fix <fix>` and :doc:`Compute <compute>` base classes.
- The Group class manipulates groups that atoms are assigned to via the
:doc:`group <group>` command. It also has functions to compute
various attributes of groups of atoms.
- The Output class is used to generate 3 kinds of output from a LAMMPS
simulation: thermodynamic information printed to the screen and log
file, dump file snapshots, and restart files. These correspond to the
:doc:`Thermo <thermo_style>`, :doc:`Dump <dump>`, and
:doc:`WriteRestart <write_restart>` classes respectively. The Dump
class is a base class with several derived classes implementing
various dump style variants.
- The Timer class logs timing information, output at the end
of a run.
.. TODO section on "Spatial decomposition and parallel operations"
.. diagram of 3d processor grid, brick vs. tiled. local vs. ghost
.. atoms, 6-way communication with pack/unpack functions,
.. PBC as part of the communication
.. TODO section on "Fixes, Computes, and Variables"
.. how and when data is computed and provided and how it is
.. referenced. flags in Fix/Compute/Variable classes tell
.. style and amount of available data.

417
doc/src/pg_dev_utils.rst Normal file
View File

@ -0,0 +1,417 @@
LAMMPS utility functions
========================
The ``utils`` sub-namespace inside the ``LAMMPS_NS`` namespace provides
a collection of convenience functions and utilities that perform common
tasks that are required repeatedly throughout the LAMMPS code like
reading or writing to files with error checking or translation of
strings into specific types of numbers with checking for validity. This
reduces redundant implementations and encourages consistent behavior.
I/O with status check
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
These are wrappers around the corresponding C library calls like
``fgets()`` or ``fread()``. They will check if there were errors
on reading or an unexpected end-of-file state was reached. In that
case, the functions will stop the calculation with an error message,
indicating the name of the problematic file, if possible.
----------
.. doxygenfunction:: sfgets
:project: progguide
.. doxygenfunction:: sfread
:project: progguide
String to number conversions with validity check
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
These functions should be used to convert strings to numbers. They are
are strongly preferred over C library calls like ``atoi()`` or
``atof()`` since they check if the **entire** provided string is a valid
(floating-point or integer) number, and will error out instead of
silently returning the result of a partial conversion or zero in cases
where the string is not a valid number. This behavior allows to more
easily detect typos or issues when processing input files.
The *do_abort* flag should be set to ``true`` in case this function
is called only on a single MPI rank, as that will then trigger the
a call to ``Error::one()`` for errors instead of ``Error::all()``
and avoids a "hanging" calculation when run in parallel.
Please also see :cpp:func:`is_integer` and :cpp:func:`is_double` for
testing strings for compliance without conversion.
----------
.. doxygenfunction:: numeric
:project: progguide
.. doxygenfunction:: inumeric
:project: progguide
.. doxygenfunction:: bnumeric
:project: progguide
.. doxygenfunction:: tnumeric
:project: progguide
String processing
^^^^^^^^^^^^^^^^^
The following are functions to help with processing strings
and parsing files or arguments.
----------
.. doxygenfunction:: trim
:project: progguide
.. doxygenfunction:: trim_comment
:project: progguide
.. doxygenfunction:: count_words(const char *text)
:project: progguide
.. doxygenfunction:: count_words(const std::string &text)
:project: progguide
.. doxygenfunction:: count_words(const std::string &text, const std::string &separators)
:project: progguide
.. doxygenfunction:: trim_and_count_words
:project: progguide
.. doxygenfunction:: split_words
:project: progguide
.. doxygenfunction:: strmatch
:project: progguide
.. doxygenfunction:: is_integer
:project: progguide
.. doxygenfunction:: is_double
:project: progguide
File and path functions
^^^^^^^^^^^^^^^^^^^^^^^^^
.. doxygenfunction:: guesspath
:project: progguide
.. doxygenfunction:: path_basename
:project: progguide
.. doxygenfunction:: path_join
:project: progguide
.. doxygenfunction:: file_is_readable
:project: progguide
Potential file functions
^^^^^^^^^^^^^^^^^^^^^^^^
.. doxygenfunction:: get_potential_file_path
:project: progguide
.. doxygenfunction:: get_potential_date
:project: progguide
.. doxygenfunction:: get_potential_units
:project: progguide
.. doxygenfunction:: get_supported_conversions
:project: progguide
.. doxygenfunction:: get_conversion_factor
:project: progguide
.. doxygenfunction:: open_potential(const std::string &name, LAMMPS *lmp, int *auto_convert)
:project: progguide
Argument processing
^^^^^^^^^^^^^^^^^^^
.. doxygenfunction:: bounds
:project: progguide
.. doxygenfunction:: expand_args
:project: progguide
Convenience functions
^^^^^^^^^^^^^^^^^^^^^
.. doxygenfunction:: logmesg
:project: progguide
.. doxygenfunction:: getsyserror
:project: progguide
.. doxygenfunction:: check_packages_for_style
:project: progguide
.. doxygenfunction:: timespec2seconds
:project: progguide
.. doxygenfunction:: date2num
:project: progguide
Customized standard functions
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
.. doxygenfunction:: merge_sort
:project: progguide
---------------------------
Tokenizer classes
=================
The purpose of the tokenizer classes is to simplify the recurring task
of breaking lines of text down into words and/or numbers.
Traditionally, LAMMPS code would be using the ``strtok()`` function from
the C library for that purpose, but that function has two significant
disadvantages: 1) it cannot be used concurrently from different LAMMPS
instances since it stores its status in a global variable and 2) it
modifies the string that it is processing. These classes were
implemented to avoid both of these issues and also to reduce the amount
of code that needs to be written.
The basic procedure is to create an instance of the tokenizer class with
the string to be processed as an argument and then do a loop until all
available tokens are read. The constructor has a default set of
separator characters, but that can be overridden. The default separators
are all "whitespace" characters, i.e. the space character, the tabulator
character, the carriage return character, the linefeed character, and
the form feed character.
.. code-block:: C++
:caption: Tokenizer class example listing entries of the PATH environment variable
#include "tokenizer.h"
#include <cstdlib>
#include <string>
#include <iostream>
using namespace LAMMPS_NS;
int main(int, char **)
{
const char *path = getenv("PATH");
if (path != nullptr) {
Tokenizer p(path,":");
while (p.has_next())
std::cout << "Entry: " << p.next() << "\n";
}
return 0;
}
Most tokenizer operations cannot fail except for
:cpp:func:`LAMMPS_NS::Tokenizer::next` (when used without first
checking with :cpp:func:`LAMMPS_NS::Tokenizer::has_next`) and
:cpp:func:`LAMMPS_NS::Tokenizer::skip`. In case of failure, the class
will throw an exception, so you may need to wrap the code using the
tokenizer into a ``try`` / ``catch`` block to handle errors. The
:cpp:class:`LAMMPS_NS::ValueTokenizer` class may also throw an exception
when a (type of) number is requested as next token that is not
compatible with the string representing the next word.
.. code-block:: C++
:caption: ValueTokenizer class example with exception handling
#include "tokenizer.h"
#include <cstdlib>
#include <string>
#include <iostream>
using namespace LAMMPS_NS;
int main(int, char **)
{
const char *text = "1 2 3 4 5 20.0 21 twentytwo 2.3";
double num1(0),num2(0),num3(0),num4(0);
ValueTokenizer t(text);
// read 4 doubles after skipping over 5 numbers
try {
t.skip(5);
num1 = t.next_double();
num2 = t.next_double();
num3 = t.next_double();
num4 = t.next_double();
} catch (TokenizerException &e) {
std::cout << "Reading numbers failed: " << e.what() << "\n";
}
std::cout << "Values: " << num1 << " " << num2 << " " << num3 << " " << num4 << "\n";
return 0;
}
This code example should produce the following output:
.. code-block::
Reading numbers failed: Not a valid floating-point number: 'twentytwo'
Values: 20 21 0 0
----------
.. doxygenclass:: LAMMPS_NS::Tokenizer
:project: progguide
:members:
.. doxygenclass:: LAMMPS_NS::TokenizerException
:project: progguide
:members:
.. doxygenclass:: LAMMPS_NS::ValueTokenizer
:project: progguide
:members:
.. doxygenclass:: LAMMPS_NS::InvalidIntegerException
:project: progguide
:members: what
.. doxygenclass:: LAMMPS_NS::InvalidFloatException
:project: progguide
:members: what
File reader classes
====================
The purpose of the file reader classes is to simplify the recurring task
of reading and parsing files. They can use the
:cpp:class:`LAMMPS_NS::ValueTokenizer` class to process the read in
text. The :cpp:class:`LAMMPS_NS::TextFileReader` is a more general
version while :cpp:class:`LAMMPS_NS::PotentialFileReader` is specialized
to implement the behavior expected for looking up and reading/parsing
files with potential parameters in LAMMPS. The potential file reader
class requires a LAMMPS instance, requires to be run on MPI rank 0 only,
will use the :cpp:func:`LAMMPS_NS::utils::get_potential_file_path`
function to look up and open the file, and will call the
:cpp:class:`LAMMPS_NS::Error` class in case of failures to read or to
convert numbers, so that LAMMPS will be aborted.
.. code-block:: C++
:caption: Use of PotentialFileReader class in pair style coul/streitz
PotentialFileReader reader(lmp, file, "coul/streitz");
char * line;
while((line = reader.next_line(NPARAMS_PER_LINE))) {
try {
ValueTokenizer values(line);
std::string iname = values.next_string();
int ielement;
for (ielement = 0; ielement < nelements; ielement++)
if (iname == elements[ielement]) break;
if (nparams == maxparam) {
maxparam += DELTA;
params = (Param *) memory->srealloc(params,maxparam*sizeof(Param),
"pair:params");
}
params[nparams].ielement = ielement;
params[nparams].chi = values.next_double();
params[nparams].eta = values.next_double();
params[nparams].gamma = values.next_double();
params[nparams].zeta = values.next_double();
params[nparams].zcore = values.next_double();
} catch (TokenizerException & e) {
error->one(FLERR, e.what());
}
nparams++;
}
A file that would be parsed by the reader code fragment looks like this:
# DATE: 2015-02-19 UNITS: metal CONTRIBUTOR: Ray Shan CITATION: Streitz and Mintmire, Phys Rev B, 50, 11996-12003 (1994)
#
# X (eV) J (eV) gamma (1/\AA) zeta (1/\AA) Z (e)
Al 0.000000 10.328655 0.000000 0.968438 0.763905
O 5.484763 14.035715 0.000000 2.143957 0.000000
----------
.. doxygenclass:: LAMMPS_NS::TextFileReader
:project: progguide
:members:
.. doxygenclass:: LAMMPS_NS::PotentialFileReader
:project: progguide
:members:
----------
Memory pool classes
===================
The memory pool classes are used for cases where otherwise many
small memory allocations would be needed and where the data would
be either all used or all freed. One example for that is the
storage of neighbor lists. The memory management strategy is
based on the assumption that allocations will be in chunks of similar
sizes. The allocation is then not done per individual call for a
reserved chunk of memory, but for a "page" that can hold multiple
chunks of data. A parameter for the maximum chunk size must be
provided, as that is used to determine whether a new page of memory
must be used.
The :cpp:class:`MyPage <LAMMPS_NS::MyPage>` class offers two ways to
reserve a chunk: 1) with :cpp:func:`get() <LAMMPS_NS::MyPage::get>` the
chunk size needs to be known in advance, 2) with :cpp:func:`vget()
<LAMMPS_NS::MyPage::vget>` a pointer to the next chunk is returned, but
its size is registered later with :cpp:func:`vgot()
<LAMMPS_NS::MyPage::vgot>`.
.. code-block:: C++
:caption: Example of using :cpp:class:`MyPage <LAMMPS_NS::MyPage>`
#include "my_page.h"
using namespace LAMMPS_NS;
MyPage<double> *dpage = new MyPage<double>;
// max size of chunk: 256, size of page: 10240 doubles (=81920 bytes)
dpage->init(256,10240);
double **build_some_lists(int num)
{
dpage->reset();
double **dlist = new double*[num];
for (int i=0; i < num; ++i) {
double *dptr = dpage.vget();
int jnum = 0;
for (int j=0; j < jmax; ++j) {
// compute some dvalue for eligible loop index j
dptr[j] = dvalue;
++jnum;
}
if (dpage.status() != 0) {
// handle out of memory or jnum too large errors
}
dpage.vgot(jnum);
dlist[i] = dptr;
}
return dlist;
}
----------
.. doxygenclass:: LAMMPS_NS::MyPage
:project: progguide
:members:
.. doxygenclass:: LAMMPS_NS::MyPoolChunk
:project: progguide
:members:

253
doc/src/pg_dev_write.rst Normal file
View File

@ -0,0 +1,253 @@
Writing LAMMPS styles
=====================
The :doc:`Modify` section of the manual gives an overview of how LAMMPS can
be extended by writing new classes that derive from existing
parent classes in LAMMPS. Here, some specific coding
details are provided for writing code for LAMMPS.
Writing a new fix style
^^^^^^^^^^^^^^^^^^^^^^^
Writing fixes is a flexible way of extending LAMMPS. Users can
implement many things using fixes:
- changing particles attributes (positions, velocities, forces, etc.). Examples: FixNVE, FixFreeze.
- reading/writing data. Example: FixRestart.
- adding or modifying properties due to geometry. Example: FixWall.
- interacting with other subsystems or external code: Examples: FixTTM, FixExternal, FixLATTE
- saving information for analysis or future use (previous positions,
for instance). Examples: Fix AveTime, FixStoreState.
All fixes are derived from the Fix base class and must have a
constructor with the signature: ``FixPrintVel(class LAMMPS *, int, char **)``.
Every fix must be registered in LAMMPS by writing the following lines
of code in the header before include guards:
.. code-block:: c
#ifdef FIX_CLASS
FixStyle(print/vel,FixPrintVel)
#else
/* the definition of the FixPrintVel class comes here */
...
#endif
Where ``print/vel`` is the style name of your fix in the input script and
``FixPrintVel`` is the name of the class. The header file would be called
``fix_print_vel.h`` and the implementation file ``fix_print_vel.cpp``.
These conventions allow LAMMPS to automatically integrate it into the
executable when compiling and associate your new fix class with the designated
keyword when it parses the input script.
Let's write a simple fix which will print the average velocity at the end
of each timestep. First of all, implement a constructor:
.. code-block:: C++
FixPrintVel::FixPrintVel(LAMMPS *lmp, int narg, char **arg)
: Fix(lmp, narg, arg)
{
if (narg < 4)
error->all(FLERR,"Illegal fix print/vel command");
nevery = force->inumeric(FLERR,arg[3]);
if (nevery <= 0)
error->all(FLERR,"Illegal fix print/vel command");
}
In the constructor you should parse your fix arguments which are
specified in the script. All fixes have pretty the same syntax:
``fix <fix-ID> <fix group> <fix name> <fix arguments ...>``. The
first 3 parameters are parsed by Fix base class constructor, while
``<fix arguments>`` should be parsed by you. In our case, we need to
specify how often we want to print an average velocity. For instance,
once in 50 timesteps: ``fix 1 print/vel 50``. There is a special variable
in the Fix class called ``nevery`` which specifies how often the method
``end_of_step()`` is called. Thus all we need to do is just set it up.
The next method we need to implement is ``setmask()``:
.. code-block:: C++
int FixPrintVel::setmask()
{
int mask = 0;
mask |= FixConst::END_OF_STEP;
return mask;
}
Here the user specifies which methods of your fix should be called
during execution. The constant ``END_OF_STEP`` corresponds to the
``end_of_step()`` method. The most important available methods that
are called during a timestep and the order in which they are called
are shown in the previous section.
.. code-block:: C++
void FixPrintVel::end_of_step()
{
// for add3, scale3
using namespace MathExtra;
double** v = atom->v;
int nlocal = atom->nlocal;
double localAvgVel[4]; // 4th element for particles count
memset(localAvgVel, 0, 4 * sizeof(double));
for (int particleInd = 0; particleInd < nlocal; ++particleInd) {
add3(localAvgVel, v[particleInd], localAvgVel);
}
localAvgVel[3] = nlocal;
double globalAvgVel[4];
memset(globalAvgVel, 0, 4 * sizeof(double));
MPI_Allreduce(localAvgVel, globalAvgVel, 4, MPI_DOUBLE, MPI_SUM, world);
scale3(1.0 / globalAvgVel[3], globalAvgVel);
if ((comm->me == 0) && screen) {
fmt::print(screen,"{}, {}, {}\n",
globalAvgVel[0], globalAvgVel[1], globalAvgVel[2]);
}
}
In the code above, we use MathExtra routines defined in
``math_extra.h``. There are bunch of math functions to work with
arrays of doubles as with math vectors. It is also important to note
that LAMMPS code should always assume to be run in parallel and that
atom data is thus distributed across the MPI ranks. Thus you can
only process data from local atoms directly and need to use MPI library
calls to combine or exchange data. For serial execution, LAMMPS
comes bundled with the MPI STUBS library that contains the MPI library
function calls in dummy versions that only work for a single MPI rank.
In this code we use an instance of Atom class. This object is stored
in the Pointers class (see ``pointers.h``) which is the base class of
the Fix base class. This object contains references to various class
instances (the original instances are created and held by the LAMMPS
class) with all global information about the simulation system.
Data from the Pointers class is available to all classes inherited from
it using protected inheritance. Hence when you write you own class,
which is going to use LAMMPS data, don't forget to inherit from Pointers
or pass an Pointer to it to all functions that need access. When writing
fixes we inherit from class Fix which is inherited from Pointers so
there is no need to inherit from it directly.
The code above computes average velocity for all particles in the
simulation. Yet you have one unused parameter in fix call from the
script: ``group_name``. This parameter specifies the group of atoms
used in the fix. So we should compute average for all particles in the
simulation only if ``group_name == "all"``, but it can be any group.
The group membership information of an atom is contained in the *mask*
property of and atom and the bit corresponding to a given group is
stored in the groupbit variable which is defined in Fix base class:
.. code-block:: C++
for (int i = 0; i < nlocal; ++i) {
if (atom->mask[i] & groupbit) {
// Do all job here
}
}
Class Atom encapsulates atoms positions, velocities, forces, etc. User
can access them using particle index. Note, that particle indexes are
usually changed every few timesteps because of neighbor list rebuilds
and spatial sorting (to improve cache efficiency).
Let us consider another Fix example: We want to have a fix which stores
atoms position from previous time step in your fix. The local atoms
indexes may not be valid on the next iteration. In order to handle
this situation there are several methods which should be implemented:
- ``double memory_usage()``: return how much memory the fix uses (optional)
- ``void grow_arrays(int)``: do reallocation of the per particle arrays in your fix
- ``void copy_arrays(int i, int j, int delflag)``: copy i-th per-particle
information to j-th. Used when atom sorting is performed. if delflag is set
and atom j owns a body, move the body information to atom i.
- ``void set_arrays(int i)``: sets i-th particle related information to zero
Note, that if your class implements these methods, it must call add calls of
add_callback and delete_callback to constructor and destructor. Since we want
to store positions of atoms from previous timestep, we need to add
``double** xold`` to the header file. Than add allocation code
to the constructor:
.. code-block:: C++
FixSavePos::FixSavePos(LAMMPS *lmp, int narg, char **arg), xold(nullptr)
{
//...
memory->create(xold, atom->nmax, 3, "FixSavePos:x");
atom->add_callback(0);
}
FixSavePos::~FixSavePos() {
atom->delete_callback(id, 0);
memory->destroy(xold);
}
Implement the aforementioned methods:
.. code-block:: C++
double FixSavePos::memory_usage()
{
int nmax = atom->nmax;
double bytes = 0.0;
bytes += nmax * 3 * sizeof(double);
return bytes;
}
void FixSavePos::grow_arrays(int nmax)
{
memory->grow(xold, nmax, 3, "FixSavePos:xold");
}
void FixSavePos::copy_arrays(int i, int j, int delflag)
{
memcpy(xold[j], xold[i], sizeof(double) * 3);
}
void FixSavePos::set_arrays(int i)
{
memset(xold[i], 0, sizeof(double) * 3);
}
int FixSavePos::pack_exchange(int i, double *buf)
{
int m = 0;
buf[m++] = xold[i][0];
buf[m++] = xold[i][1];
buf[m++] = xold[i][2];
return m;
}
int FixSavePos::unpack_exchange(int nlocal, double *buf)
{
int m = 0;
xold[nlocal][0] = buf[m++];
xold[nlocal][1] = buf[m++];
xold[nlocal][2] = buf[m++];
return m;
}
Now, a little bit about memory allocation. We use the Memory class which
is just a bunch of template functions for allocating 1D and 2D
arrays. So you need to add include ``memory.h`` to have access to them.
Finally, if you need to write/read some global information used in
your fix to the restart file, you might do it by setting flag
``restart_global = 1`` in the constructor and implementing methods void
``write_restart(FILE *fp)`` and ``void restart(char *buf)``.
If, in addition, you want to write the per-atom property to restart
files additional settings and functions are needed:
- a fix flag indicating this needs to be set ``restart_peratom = 1;``
- ``atom->add_callback()`` and ``atom->delete_callback()`` must be called
a second time with the final argument set to 1 instead of 0 (indicating
restart processing instead of per-atom data memory management).
- the functions ``void pack_restart(int i, double *buf)`` and
``void unpack_restart(int nlocal, int nth)`` need to be implemented

File diff suppressed because it is too large Load Diff

7
doc/src/pg_input.rst Normal file
View File

@ -0,0 +1,7 @@
LAMMPS Input Base Class
************************
.. doxygenclass:: LAMMPS_NS::Input
:project: progguide
:members:

22
doc/src/pg_lammps.rst Normal file
View File

@ -0,0 +1,22 @@
LAMMPS Class
************
The LAMMPS class is encapsulating an MD simulation state and thus it is
the class that needs to be created when starting a new simulation system
state. The LAMMPS executable essentially creates one instance of this
class and passes the command line flags and tells it to process the
provided input (a file or ``stdin``). It shuts the class down when
control is returned to it and then exits. When using LAMMPS as a
library from another code it is required to create an instance of this
class, either directly from C++ with ``new LAMMPS()`` or through one
of the library interface functions like :cpp:func:`lammps_open` of the
C-library interface, or the :py:class:`lammps.lammps` class constructor
of the Python module, or the :f:func:`lammps` constructor of the Fortran
module.
--------------------
.. doxygenclass:: LAMMPS_NS::LAMMPS
:project: progguide
:members:

View File

@ -3,3 +3,7 @@
display: block;
margin-bottom: 0.809em;
}
.versionmodified {
font-weight: bold;
}

View File

@ -1300,6 +1300,7 @@ initializations
initio
InP
inregion
instantiation
Institut
integrators
Integrators
@ -1752,6 +1753,7 @@ Mattox
Mattson
maxangle
maxbond
maxchunk
maxelt
maxeval
maxfiles
@ -1975,6 +1977,7 @@ MxN
myCompute
myIndex
mylammps
MyPool
mysocket
myTemp
myVec

View File

@ -54,12 +54,11 @@ MODULE LIBLAMMPS
! interface definitions for calling functions in library.cpp
INTERFACE
FUNCTION lammps_open(argc,argv,comm,handle) &
FUNCTION lammps_open(argc,argv,comm) &
BIND(C, name='lammps_open_fortran')
IMPORT :: c_ptr, c_int
INTEGER(c_int), VALUE, INTENT(in) :: argc, comm
TYPE(c_ptr), DIMENSION(*), INTENT(in) :: argv
TYPE(c_ptr), INTENT(out) :: handle
TYPE(c_ptr) :: lammps_open
END FUNCTION lammps_open
@ -161,7 +160,7 @@ CONTAINS
ENDIF
IF (PRESENT(comm)) THEN
lmp_open%handle = lammps_open(argc,argv,comm,dummy)
lmp_open%handle = lammps_open(argc,argv,comm)
ELSE
lmp_open%handle = lammps_open_no_mpi(argc,argv,dummy)
END IF

View File

@ -49,6 +49,8 @@ packages_ntopo.h
# other auto-generated files
lmpinstalledpkgs.h
lmpgitversion.h
# removed on 9 Sep 2020
mergesort.h
# renamed on 8 May 2020
fix_meso.cpp
fix_meso.h

View File

@ -164,7 +164,7 @@ FixHyperLocal::FixHyperLocal(LAMMPS *lmp, int narg, char **arg) :
maxbondperatom = FCCBONDS;
numcoeff = NULL;
clist = NULL;
cpage = new MyPage<OneCoeff>;
cpage = new MyPage<HyperOneCoeff>;
cpage->init(maxbondperatom,1024*maxbondperatom,1);
// set comm sizes needed by this fix
@ -976,7 +976,7 @@ void FixHyperLocal::build_bond_list(int natom)
memory->sfree(clist);
maxcoeff = atom->nmax;
memory->create(numcoeff,maxcoeff,"hyper/local:numcoeff");
clist = (OneCoeff **) memory->smalloc(maxcoeff*sizeof(OneCoeff *),
clist = (HyperOneCoeff **) memory->smalloc(maxcoeff*sizeof(HyperOneCoeff *),
"hyper/local:clist");
}
@ -1741,7 +1741,7 @@ double FixHyperLocal::memory_usage()
bytes += 2*maxall * sizeof(double); // maxstrain,maxstrain_domain
if (checkbias) bytes += maxall * sizeof(tagint); // biasflag
bytes += maxcoeff * sizeof(int); // numcoeff
bytes += maxcoeff * sizeof(OneCoeff *); // clist
bytes += maxlocal*maxbondperatom * sizeof(OneCoeff); // cpage estimate
bytes += maxcoeff * sizeof(HyperOneCoeff *); // clist
bytes += maxlocal*maxbondperatom * sizeof(HyperOneCoeff); // cpage estimate
return bytes;
}

View File

@ -23,6 +23,8 @@ FixStyle(hyper/local,FixHyperLocal)
#include "fix_hyper.h"
namespace LAMMPS_NS {
// forward declaration. struct HyperOneCoeff is defined in my_page.h
struct HyperOneCoeff;
class FixHyperLocal : public FixHyper {
public:
@ -183,13 +185,8 @@ class FixHyperLocal : public FixHyper {
// data structs for persisting bias coeffs when bond list is reformed
struct OneCoeff {
double biascoeff;
tagint tag;
};
MyPage<OneCoeff> *cpage; // pages of OneCoeff datums for clist
OneCoeff **clist; // ptrs to vectors of bias coeffs for each atom
MyPage<HyperOneCoeff> *cpage;// pages of OneCoeff datums for clist
HyperOneCoeff **clist; // ptrs to vectors of bias coeffs for each atom
int *numcoeff; // # of bias coeffs per atom (one per bond)
int maxcoeff; // allocate sized of clist and numcoeff

View File

@ -47,6 +47,33 @@ using namespace MathConst;
/* ---------------------------------------------------------------------- */
/** \class LAMMPS_NS::Atom
* \brief Class to provide access to atom data
\verbatim embed:rst
The Atom class provides access to atom style related global settings and
per-atom data that is stored with atoms and migrates with them from
sub-domain to sub-domain as atoms move around. This includes topology
data, which is stored with either one specific atom or all atoms involved
depending on the settings of the :doc:`newton command <newton>`.
The actual per-atom data is allocated and managed by one of the various
classes derived from the AtomVec class as determined by
the :doc:`atom_style command <atom_style>`. The pointers in the Atom class
are updated by the AtomVec class as needed.
\endverbatim
*/
/** Atom class constructor
*
* This resets and initialized all kinds of settings,
* parameters, and pointer variables for per-atom arrays.
* This also initializes the factory for creating
* instances of classes derived from the AtomVec base
* class, which correspond to the selected atom style.
*
* \param lmp pointer to the base LAMMPS class */
Atom::Atom(LAMMPS *lmp) : Pointers(lmp)
{
natoms = 0;
@ -688,7 +715,6 @@ AtomVec *Atom::avec_creator(LAMMPS *lmp)
return new T(lmp);
}
/* ---------------------------------------------------------------------- */
void Atom::init()
@ -2301,12 +2327,17 @@ int Atom::find_custom(const char *name, int &flag)
return -1;
}
/* ----------------------------------------------------------------------
add a custom variable with name of type flag = 0/1 for int/double
assumes name does not already exist
return index in ivector or dvector of its location
------------------------------------------------------------------------- */
/** \brief Add a custom per-atom property with the given name and type
\verbatim embed:rst
This function will add a custom per-atom property with the name "name"
as either list of int or double to the list of custom properties. This
function is called, e.g. from :doc:`fix property/atom <fix_property_atom>`.
\endverbatim
* \param name Name of the property (w/o a "d_" or "i_" prefix)
* \param flag Data type of property: 0 for int, 1 for double
* \return Index of property in the respective list of properties
*/
int Atom::add_custom(const char *name, int flag)
{
int index;
@ -2338,12 +2369,19 @@ int Atom::add_custom(const char *name, int flag)
return index;
}
/* ----------------------------------------------------------------------
remove a custom variable of type flag = 0/1 for int/double at index
free memory for vector and name and set ptrs to NULL
ivector/dvector and iname/dname lists never shrink
------------------------------------------------------------------------- */
/*! \brief Remove a custom per-atom property of a given type
*
\verbatim embed:rst
This will remove a property that was requested e.g. by the
:doc:`fix property/atom <fix_property_atom>` command. It frees the
allocated memory and sets the pointer to ``NULL`` to the entry in
the list can be reused. The lists of those pointers will never be
compacted or never shrink, so that index to name mappings remain valid.
\endverbatim
*
* \param flag whether the property is integer (=0) or double (=1)
* \param index of that property in the respective list.
*/
void Atom::remove_custom(int flag, int index)
{
if (flag == 0) {
@ -2359,16 +2397,123 @@ void Atom::remove_custom(int flag, int index)
}
}
/* ----------------------------------------------------------------------
return a pointer to a named internal variable
if don't recognize name, return NULL
------------------------------------------------------------------------- */
/** Provide access to internal data of the Atom class by keyword
*
\verbatim embed:rst
This function is a way to access internal per-atom data. This data is
distributed across MPI ranks and thus only the data for "local" atoms
are expected to be available. Whether also data for "ghost" atoms is
stored and up-to-date depends on various simulation settings.
This table lists a large part of the supported names, their data types,
length of the data area, and a short description.
.. list-table::
:header-rows: 1
:widths: auto
* - Name
- Type
- Items per atom
- Description
* - mass
- double
- 1
- per-type mass. This array is **NOT** a per-atom array
but of length ``ntypes+1``, element 0 is ignored.
* - id
- tagint
- 1
- atom ID of the particles
* - type
- int
- 1
- atom type of the particles
* - mask
- int
- 1
- bitmask for mapping to groups. Individual bits are set
to 0 or 1 for each group.
* - image
- imageint
- 1
- 3 image flags encoded into a single integer.
See :cpp:func:`lammps_encode_image_flags`.
* - x
- double
- 3
- x-, y-, and z-coordinate of the particles
* - v
- double
- 3
- x-, y-, and z-component of the velocity of the particles
* - f
- double
- 3
- x-, y-, and z-component of the force on the particles
* - molecule
- int
- 1
- molecule ID of the particles
* - q
- double
- 1
- charge of the particles
* - mu
- double
- 3
- dipole moment of the particles
* - omega
- double
- 3
- x-, y-, and z-component of rotational velocity of the particles
* - angmom
- double
- 3
- x-, y-, and z-component of angular momentum of the particles
* - torque
- double
- 3
- x-, y-, and z-component of the torque on the particles
* - radius
- double
- 1
- radius of the (extended) particles
* - rmass
- double
- 1
- per-atom mass of the particles. ``NULL`` if per-type masses are
used. See the :cpp:func:`rmass_flag<lammps_extract_setting>` setting.
* - ellipsoid
- int
- 1
- 1 if the particle is an ellipsoidal particle, 0 if not
* - line
- int
- 1
- 1 if the particle is a line particle, 0 if not
* - tri
- int
- 1
- 1 if the particle is a triangulated particle, 0 if not
* - body
- int
- 1
- 1 if the particle is a body particle, 0 if not
\endverbatim
*
* \param name string with the keyword of the desired property.
Typically the name of the pointer variable returned
* \return pointer to the requested data cast to ``void *`` or NULL */
void *Atom::extract(char *name)
{
// --------------------------------------------------------------------
// 4th customization section: customize by adding new variable name
/* NOTE: this array is only of length ntypes+1 */
if (strcmp(name,"mass") == 0) return (void *) mass;
if (strcmp(name,"id") == 0) return (void *) tag;
@ -2389,6 +2534,7 @@ void *Atom::extract(char *name)
if (strcmp(name,"ellipsoid") == 0) return (void *) ellipsoid;
if (strcmp(name,"line") == 0) return (void *) line;
if (strcmp(name,"tri") == 0) return (void *) tri;
if (strcmp(name,"body") == 0) return (void *) body;
if (strcmp(name,"vfrac") == 0) return (void *) vfrac;
if (strcmp(name,"s0") == 0) return (void *) s0;

View File

@ -555,7 +555,7 @@ bigint AtomVecBody::memory_usage_bonus()
{
bigint bytes = 0;
bytes += nmax_bonus*sizeof(Bonus);
bytes += icp->size + dcp->size;
bytes += icp->size() + dcp->size();
int nall = nlocal_bonus + nghost_bonus;
for (int i = 0; i < nall; i++) {

View File

@ -32,8 +32,6 @@ using namespace LAMMPS_NS;
#if defined(LMP_QSORT)
// allocate space for static class variable
Dump *Dump::dumpptr;
#else
#include "mergesort.h"
#endif
#define BIG 1.0e20
@ -766,9 +764,9 @@ void Dump::sort()
#else
if (!reorderflag) {
for (i = 0; i < nme; i++) index[i] = i;
if (sortcol == 0) merge_sort(index,nme,(void *)this,idcompare);
else if (sortorder == ASCEND) merge_sort(index,nme,(void *)this,bufcompare);
else merge_sort(index,nme,(void *)this,bufcompare_reverse);
if (sortcol == 0) utils::merge_sort(index,nme,(void *)this,idcompare);
else if (sortorder == ASCEND) utils::merge_sort(index,nme,(void *)this,bufcompare);
else utils::merge_sort(index,nme,(void *)this,bufcompare_reverse);
}
#endif

View File

@ -60,6 +60,37 @@ using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */
/** \class LAMMPS_NS::Input
* \brief Class for processing commands and input files
*
\verbatim embed:rst
The Input class contains methods for reading, pre-processing and
parsing LAMMPS commands and input files and will dispatch commands
to the respective class instances or contains the code to execute
the commands directly. It also contains the instance of the
Variable class which performs computations and text substitutions.
\endverbatim */
/** Input class constructor
*
\verbatim embed:rst
This sets up the input processing, processes the *-var* and *-echo*
command line flags, holds the factory of commands and creates and
initializes an instance of the Variable class.
To execute a command, a specific class instance, derived from
:cpp:class:`Pointers`, is created, then its ``command()`` member
function executed, and finally the class instance is deleted.
\endverbatim
*
* \param lmp pointer to the base LAMMPS class
* \param argc number of entries in *argv*
* \param argv argument vector */
Input::Input(LAMMPS *lmp, int argc, char **argv) : Pointers(lmp)
{
MPI_Comm_rank(world,&me);
@ -137,10 +168,15 @@ Input::~Input()
delete command_map;
}
/* ----------------------------------------------------------------------
process all input from infile
infile = stdin or file if command-line arg "-in" was used
------------------------------------------------------------------------- */
/** Process all input from the ``FILE *`` pointer *infile*
*
\verbatim embed:rst
This will read lines from *infile*, parse and execute them until the end
of the file is reached. The *infile* pointer will usually point to
``stdin`` or the input file given with the ``-in`` command line flag.
\endverbatim */
void Input::file()
{
@ -229,10 +265,21 @@ void Input::file()
}
}
/* ----------------------------------------------------------------------
process all input from file at filename
mostly called from library interface
------------------------------------------------------------------------- */
/** Process all input from the file *filename*
*
\verbatim embed:rst
This function opens the file at the path *filename*, put the current
file pointer stored in *infile* on a stack and instead assign *infile*
with the newly opened file pointer. Then it will call the
:cpp:func:`Input::file() <LAMMPS_NS::Input::file()>` function to read,
parse and execute the contents of that file. When the end of the file
is reached, it is closed and the previous file pointer from the infile
file pointer stack restored to *infile*.
\endverbatim
*
* \param filename name of file with LAMMPS commands */
void Input::file(const char *filename)
{
@ -263,11 +310,19 @@ void Input::file(const char *filename)
}
}
/* ----------------------------------------------------------------------
invoke one command in single
first copy to line, then parse, then execute it
return command name to caller
------------------------------------------------------------------------- */
/** Process a single command from a string in *single*
*
\verbatim embed:rst
This function takes the text in *single*, makes a copy, parses that,
executes the command and returns the name of the command (without the
arguments). If there was no command in *single* it will return
``NULL``.
\endverbatim
*
* \param single string with LAMMPS command
* \return string with name of the parsed command w/o arguments */
char *Input::one(const std::string &single)
{

View File

@ -31,7 +31,6 @@ using namespace LAMMPS_NS;
int *Irregular::proc_recv_copy;
static int compare_standalone(const void *, const void *);
#else
#include "mergesort.h"
// prototype for non-class function
static int compare_standalone(const int, const int, void *);
#endif
@ -441,7 +440,7 @@ int Irregular::create_atom(int n, int *sizes, int *proclist, int sortflag)
proc_recv_copy = proc_recv;
qsort(order,nrecv_proc,sizeof(int),compare_standalone);
#else
merge_sort(order,nrecv_proc,(void *)proc_recv,compare_standalone);
utils::merge_sort(order,nrecv_proc,(void *)proc_recv,compare_standalone);
#endif
int j;
@ -715,7 +714,7 @@ int Irregular::create_data(int n, int *proclist, int sortflag)
proc_recv_copy = proc_recv;
qsort(order,nrecv_proc,sizeof(int),compare_standalone);
#else
merge_sort(order,nrecv_proc,(void *)proc_recv,compare_standalone);
utils::merge_sort(order,nrecv_proc,(void *)proc_recv,compare_standalone);
#endif
int j;
@ -889,7 +888,7 @@ int Irregular::create_data_grouped(int n, int *procs, int sortflag)
proc_recv_copy = proc_recv;
qsort(order,nrecv_proc,sizeof(int),compare_standalone);
#else
merge_sort(order,nrecv_proc,(void *)proc_recv,compare_standalone);
utils::merge_sort(order,nrecv_proc,(void *)proc_recv,compare_standalone);
#endif
int j;

View File

@ -80,14 +80,30 @@ struct LAMMPS_NS::package_styles_lists {
using namespace LAMMPS_NS;
/* ----------------------------------------------------------------------
start up LAMMPS
allocate fundamental classes (memory, error, universe, input)
parse input switches
initialize communicators, screen & logfile output
input is allocated at end after MPI info is setup
------------------------------------------------------------------------- */
/** \class LAMMPS_NS::LAMMPS
* \brief LAMMPS simulation instance
*
* The LAMMPS class contains pointers of all constituent class instances
* and global variables that are used by a LAMMPS simulation. Its contents
* represent the entire state of the simulation.
*
* The LAMMPS class manages the components of an MD simulation by creating,
* deleting, and initializing instances of the classes it is composed of,
* processing command line flags, and providing access to some global properties.
* The specifics of setting up and running a simulation are handled by the
* individual component class instances. */
/** Create a LAMMPS simulation instance
*
* The LAMMPS constructor starts up a simulation by allocating all
* fundamental classes in the necessary order, parses input switches
* and their arguments, initializes communicators, screen and logfile
* output FILE pointers.
*
* \param narg number of arguments
* \param arg list of arguments
* \param communicator MPI communicator used by this LAMMPS instance
*/
LAMMPS::LAMMPS(int narg, char **arg, MPI_Comm communicator) :
memory(NULL), error(NULL), universe(NULL), input(NULL), atom(NULL),
update(NULL), neighbor(NULL), comm(NULL), domain(NULL), force(NULL),
@ -636,14 +652,13 @@ LAMMPS::LAMMPS(int narg, char **arg, MPI_Comm communicator) :
}
}
/* ----------------------------------------------------------------------
shutdown LAMMPS
delete top-level classes
close screen and log files in world and universe
output files were already closed in destroy()
delete fundamental classes
------------------------------------------------------------------------- */
/** Shut down a LAMMPS simulation instance
*
* The LAMMPS destructor shuts down the simulation by deleting top-level class
* instances, closing screen and log files for the global instance (aka "world")
* and files and MPI communicators in sub-partitions ("universes"). Then it
* deletes the fundamental class instances and copies of data inside the class.
*/
LAMMPS::~LAMMPS()
{
const int me = comm->me;
@ -989,6 +1004,11 @@ void _noopt LAMMPS::init_pkg_lists()
#undef REGION_CLASS
}
/** Return true if a LAMMPS package is enabled in this binary
*
* \param pkg name of package
* \return true if yes, else false
*/
bool LAMMPS::is_installed_pkg(const char *pkg)
{
for (int i=0; installed_packages[i] != NULL; ++i)
@ -1005,6 +1025,16 @@ bool LAMMPS::is_installed_pkg(const char *pkg)
} \
}
/** \brief Return name of package that a specific style belongs to
*
* This function checks the given name against all list of styles
* for all type of styles and if the name and the style match, it
* returns which package this style belongs to.
*
* \param style Type of style (e.g. atom, pair, fix, etc.)
* \param name Name of style
* \return Name of the package this style is part of
*/
const char *LAMMPS::match_style(const char *style, const char *name)
{
check_for_match(angle,style,name);

View File

@ -85,8 +85,10 @@ class LAMMPS {
struct package_styles_lists *pkg_lists;
void init_pkg_lists();
void help();
LAMMPS() {}; // prohibit using the default constructor
LAMMPS(const LAMMPS &) {}; // prohibit using the copy constructor
/// Default constructor. Declared private to prohibit its use
LAMMPS() {};
/// Copy constructor. Declared private to prohibit its use
LAMMPS(const LAMMPS &) {};
};
}

View File

@ -108,18 +108,23 @@ thus is otherwise ignored. However ``argc`` may be set to 0 and then
``argv`` may be ``NULL``. If MPI is not yet initialized, ``MPI_Init()``
will be called during creation of the LAMMPS class instance.
The function returns a pointer to the created LAMMPS class. If for some
reason the initialization of the LAMMPS instance fails, the function
returns ``NULL``. For backward compatibility it is also possible to
provide the address of a pointer variable as argument *ptr*\ . This
argument may be ``NULL`` and is then ignored.
If for some reason the creation or initialization of the LAMMPS instance
fails a null pointer is returned.
.. versionchanged:: 15Sep2020
This function now has the pointer to the created LAMMPS class
instance as return value. For backward compatibility it is still
possible to provide the address of a pointer variable as final
argument *ptr*\ . This use is deprecated and may be removed in
the future. The *ptr* argument may be ``NULL`` and is then ignored.
.. note::
This function is not declared when the code linking to the LAMMPS
library interface is compiled with ``-DLAMMPS_LIB_NO_MPI``, or
contains a ``#define LAMMPS_LIB_NO_MPI 1`` statement before
``#include "library.h"``. In that case, you need to use the
``#include "library.h"``. In that case, you must use the
:cpp:func:`lammps_open_no_mpi` function.
\endverbatim
@ -169,6 +174,17 @@ library was compiled in serial mode, but the calling code runs in
parallel and the ``MPI_Comm`` data type of the STUBS library would not
be compatible with that of the calling code.
If for some reason the creation or initialization of the LAMMPS instance
fails a null pointer is returned.
.. versionchanged:: 15Sep2020
This function now has the pointer to the created LAMMPS class
instance as return value. For backward compatibility it is still
possible to provide the address of a pointer variable as final
argument *ptr*\ . This use is deprecated and may be removed in
the future. The *ptr* argument may be ``NULL`` and is then ignored.
\endverbatim
*
* \param argc number of command line arguments
@ -195,20 +211,23 @@ module. Internally it converts the *f_comm* argument into a C-style MPI
communicator with ``MPI_Comm_f2c()`` and then calls
:cpp:func:`lammps_open`.
If for some reason the creation or initialization of the LAMMPS instance
fails a null pointer is returned.
.. versionadded:: 15Sep2020
\endverbatim
*
* \param argc number of command line arguments
* \param argv list of command line argument strings
* \param f_comm Fortran style MPI communicator for this LAMMPS instance
* \param ptr pointer to a void pointer variable
* which serves as a handle; may be ``NULL``
* \return pointer to new LAMMPS instance cast to ``void *`` */
void *lammps_open_fortran(int argc, char **argv, int f_comm, void **ptr)
void *lammps_open_fortran(int argc, char **argv, int f_comm)
{
lammps_mpi_init();
MPI_Comm c_comm = MPI_Comm_f2c((MPI_Fint)f_comm);
return lammps_open(argc, argv, c_comm, ptr);
return lammps_open(argc, argv, c_comm, nullptr);
}
/* ---------------------------------------------------------------------- */
@ -244,6 +263,8 @@ The MPI standard requires that any MPI application must call
calls. This function checks, whether MPI is already initialized and
calls ``MPI_Init()`` in case it is not.
.. versionadded:: 15Sep2020
\endverbatim */
void lammps_mpi_init()
@ -274,6 +295,8 @@ before exiting the program to wait until all (parallel) tasks are
completed and then MPI is cleanly shut down. After this function no
more MPI calls may be made.
.. versionadded:: 15Sep2020
\endverbatim */
void lammps_mpi_finalize()

View File

@ -77,7 +77,7 @@ extern "C" {
void *lammps_open(int argc, char **argv, MPI_Comm comm, void **ptr);
#endif
void *lammps_open_no_mpi(int argc, char **argv, void **ptr);
void *lammps_open_fortran(int argc, char **argv, int f_comm, void **ptr);
void *lammps_open_fortran(int argc, char **argv, int f_comm);
void lammps_close(void *handle);
void lammps_mpi_init();
void lammps_mpi_finalize();

View File

@ -1,124 +0,0 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#ifndef LMP_MERGESORT
#define LMP_MERGESORT
#include <cstring>
// custom hybrid upward merge sort implementation with support to pass
// an opaque pointer to the comparison function, e.g. for access to
// class members. this avoids having to use global variables.
// for improved performance, we employ an in-place insertion sort on
// chunks of up to 64 elements and switch to merge sort from then on.
// part 1. insertion sort for pre-sorting of small chunks
static void insertion_sort(int *index, int num, void *ptr,
int (*comp)(int, int, void*))
{
if (num < 2) return;
for (int i=1; i < num; ++i) {
int tmp = index[i];
for (int j=i-1; j >= 0; --j) {
if ((*comp)(index[j],tmp,ptr) > 0) {
index[j+1] = index[j];
} else {
index[j+1] = tmp;
break;
}
if (j == 0) index[0] = tmp;
}
}
}
// part 2. merge two sublists
static void do_merge(int *idx, int *buf, int llo, int lhi, int rlo, int rhi,
void *ptr, int (*comp)(int, int, void *))
{
int i = llo;
int l = llo;
int r = rlo;
while ((l < lhi) && (r < rhi)) {
if ((*comp)(buf[l],buf[r],ptr) < 0)
idx[i++] = buf[l++];
else idx[i++] = buf[r++];
}
while (l < lhi) idx[i++] = buf[l++];
while (r < rhi) idx[i++] = buf[r++];
}
// part 3: loop over sublists doubling in size with each iteration.
// pre-sort sublists with insertion sort for better performance.
static void merge_sort(int *index, int num, void *ptr,
int (*comp)(int, int, void *))
{
if (num < 2) return;
int chunk,i,j;
// do insertion sort on chunks of up to 64 elements
chunk = 64;
for (i=0; i < num; i += chunk) {
j = (i+chunk > num) ? num-i : chunk;
insertion_sort(index+i,j,ptr,comp);
}
// already done?
if (chunk >= num) return;
// continue with merge sort on the pre-sorted chunks.
// we need an extra buffer for temporary storage and two
// pointers to operate on, so we can swap the pointers
// rather than copying to the hold buffer in each pass
int *buf = new int[num];
int *dest = index;
int *hold = buf;
while (chunk < num) {
int m;
// swap hold and destination buffer
int *tmp = dest; dest = hold; hold = tmp;
// merge from hold array to destination array
for (i=0; i < num-1; i += 2*chunk) {
j = i + 2*chunk;
if (j > num) j=num;
m = i+chunk;
if (m > num) m=num;
do_merge(dest,hold,i,m,m,j,ptr,comp);
}
// copy all indices not handled by the chunked merge sort loop
for ( ; i < num ; i++ ) dest[i] = hold[i];
chunk *= 2;
}
// if the final sorted data is in buf, copy back to index
if (dest == buf) memcpy(index,buf,sizeof(int)*num);
delete[] buf;
}
#endif

193
src/my_page.cpp Normal file
View File

@ -0,0 +1,193 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#include "my_page.h"
#include <cstdlib>
#if defined(LMP_USER_INTEL) && !defined(LAMMPS_MEMALIGN) && !defined(_WIN32)
#define LAMMPS_MEMALIGN 64
#endif
using namespace LAMMPS_NS;
/** \class LAMMPS_NS::MyPage
* \brief Templated class for storing chunks of datums in pages.
*
* The size of the chunk may vary from call to call, but must be
* less or equal than the *maxchunk* setting.
* The chunks are not returnable like with malloc() (i.e. you cannot
* call free() on them individually). One can only reset and start over.
* The purpose of this class is to replace many small memory allocations
* via malloc() with a few large ones. Since the pages are never freed
* until the class is re-initialized, they can be re-used without having
* to re-allocate them by calling the reset() method.
*
* The settings *maxchunk*, *pagesize*, and *pagedelta* control
* the memory allocation strategy. The *maxchunk* value represents
* the expected largest number of items per chunk. If there is
* less space left on the current page, a new page is allocated
* for the next chunk. The *pagesize* value represents how many
* items can fit on a single page. It should have space for multiple
* chunks of size *maxchunk*. The combination of these two
* parameters determines how much memory is wasted by either switching
* to the next page too soon or allocating too large pages that never
* get properly used. It is an error, if a requested chunk is larger
* than *maxchunk*. The *pagedelta* parameter determines how many
* pages are allocated in one go. In combination with the *pagesize*
* setting, this determines how often blocks of memory get allocated
* (fewer allocations will result in faster execution).
*
* \note
* This is a template class with explicit instantiation. If the class
* is used with a new data type a new explicit instantiation may need to
* be added at the end of the file ``src/my_page.cpp`` to avoid symbol
* lookup errors. */
/** Create a class instance
*
* Need to call init() before use to define allocation settings */
template <class T>
MyPage<T>::MyPage() : ndatum(0), nchunk(0), pages(nullptr), page(nullptr),
npage(0), ipage(-1), index(-1), maxchunk(-1),
pagesize(-1), pagedelta(1), errorflag(0) {};
template <class T>
MyPage<T>::~MyPage() {
deallocate();
}
/** (Re-)initialize the set of pages and allocation parameters.
*
* This also frees all previously allocated storage and allocates
* the first page(s).
*
* \param user_maxchunk Expected maximum number of items for one chunk
* \param user_pagesize Number of items on a single memory page
* \param user_pagedelta Number of pages to allocate with one malloc
* \return 1 if there were invalid parameters, 2 if there was an allocation error or 0 if successful */
template<class T>
int MyPage<T>::init(int user_maxchunk, int user_pagesize,
int user_pagedelta) {
maxchunk = user_maxchunk;
pagesize = user_pagesize;
pagedelta = user_pagedelta;
if (maxchunk <= 0 || pagesize <= 0 || pagedelta <= 0) return 1;
if (maxchunk > pagesize) return 1;
// free storage if re-initialized
deallocate();
// initial page allocation
allocate();
if (errorflag) return 2;
reset();
return 0;
}
/** Pointer to location that can store N items.
*
* This will allocate more pages as needed.
* If the parameter *N* is larger than the *maxchunk*
* setting an error is flagged.
*
* \param n number of items for which storage is requested
* \return memory location or null pointer, if error or allocation failed */
template <class T>
T *MyPage<T>::get(int n) {
if (n > maxchunk) {
errorflag = 1;
return NULL;
}
ndatum += n;
nchunk++;
// return pointer from current page
if (index+n <= pagesize) {
int start = index;
index += n;
return &page[start];
}
// allocate new page
ipage++;
if (ipage == npage) {
allocate();
if (errorflag) return NULL;
}
page = pages[ipage];
index = n;
return &page[0];
}
/** Reset state of memory pool without freeing any memory */
template <class T>
void MyPage<T>::reset() {
ndatum = nchunk = 0;
index = ipage = 0;
page = (pages != nullptr) ? pages[ipage] : nullptr;
errorflag = 0;
}
/* ---------------------------------------------------------------------- */
template <class T>
void MyPage<T>::allocate() {
npage += pagedelta;
pages = (T **) realloc(pages,npage*sizeof(T *));
if (!pages) {
errorflag = 2;
return;
}
for (int i = npage-pagedelta; i < npage; i++) {
#if defined(LAMMPS_MEMALIGN)
void *ptr;
if (posix_memalign(&ptr, LAMMPS_MEMALIGN, pagesize*sizeof(T)))
errorflag = 2;
pages[i] = (T *) ptr;
#else
pages[i] = (T *) malloc(pagesize*sizeof(T));
if (!pages[i]) errorflag = 2;
#endif
}
}
/** Free all allocated pages of this class instance */
template <class T>
void MyPage<T>::deallocate() {
reset();
for (int i = 0; i < npage; i++) free(pages[i]);
free(pages);
pages = nullptr;
npage = 0;
}
// explicit instantiations
namespace LAMMPS_NS {
template class MyPage<int>;
template class MyPage<long>;
template class MyPage<long long>;
template class MyPage<double>;
template class MyPage<HyperOneCoeff>;
}

View File

@ -12,144 +12,41 @@
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
MyPage = templated class for storing chunks of datums in pages
chunks are not returnable, can only reset and start over
replaces many small mallocs with a few large mallocs
pages are never freed, so can reuse w/out reallocs
usage:
request one datum at a time, repeat, clear
request chunks of datums in each get() or vget(), repeat, clear
chunk size can vary from request to request
chunk size can be known in advance or registered after usage via vgot()
inputs:
template T = one datum, e.g. int, double, struct, int[3]
for int[3], access datum as ivec[i][2]
methods:
T *get() = return ptr to one datum
T *get(N) = return ptr to N datums, N < maxchunk required
T *vget() = return ptr to maxchunk datums, use as needed, then call vgot()
all gets return NULL if error encountered
vgot(N) = used N datums of previous vget(), N < maxchunk required
void init(maxchunk, pagesize, pagedelta)
define allocation params and allocate first page(s)
call right after constructor
can call again to reset allocation params and free previous pages
maxchunk = max # of datums in one chunk, default = 1
pagesize = # of datums in one page, default = 1024
should be big enough to store multiple chunks
pagedelta = # of pages to allocate at a time, default = 1
return 1 if bad params
void reset() = clear pages w/out freeing
int size() = return total size of allocated pages in bytes
int status() = return error status
0 = ok, 1 = chunksize > maxchunk, 2 = allocation error
templated class for storing chunks of datums in pages
------------------------------------------------------------------------- */
#ifndef LAMMPS_MY_PAGE_H
#define LAMMPS_MY_PAGE_H
#if defined(LMP_USER_INTEL) && !defined(LAMMPS_MEMALIGN) && !defined(_WIN32)
#define LAMMPS_MEMALIGN 64
#endif
#include "lmptype.h"
#include <cstdlib>
namespace LAMMPS_NS {
struct HyperOneCoeff {
double biascoeff;
tagint tag;
};
template<class T>
class MyPage {
public:
int ndatum; // total # of stored datums
int nchunk; // total # of stored chunks
MyPage();
virtual ~MyPage();
MyPage() {
ndatum = nchunk = 0;
pages = NULL;
npage = 0;
errorflag = 0;
}
int init(int user_maxchunk=1, int user_pagesize=1024,
int user_pagedelta=1);
// (re)initialize allocation params
// also allocate first page(s)
T *get(int n=1);
int init(int user_maxchunk = 1, int user_pagesize = 1024,
int user_pagedelta = 1) {
maxchunk = user_maxchunk;
pagesize = user_pagesize;
pagedelta = user_pagedelta;
if (maxchunk <= 0 || pagesize <= 0 || pagedelta <= 0) return 1;
if (maxchunk > pagesize) return 1;
// free any previously allocated pages
for (int i = 0; i < npage; i++) free(pages[i]);
free(pages);
// initial page allocation
ndatum = nchunk = 0;
pages = NULL;
npage = 0;
allocate();
if (errorflag) return 2;
ipage = index = 0;
page = pages[ipage];
return 0;
}
// free all allocated pages
~MyPage() {
for (int i = 0; i < npage; i++) free(pages[i]);
free(pages);
}
// get ptr to one datum
// return NULL if run out of memory
T *get() {
ndatum++;
nchunk++;
if (index < pagesize) return &page[index++];
ipage++;
if (ipage == npage) {
allocate();
if (errorflag) return NULL;
}
page = pages[ipage];
index = 0;
return &page[index++];
}
// get ptr to location that can store N datums
// error if N > maxchunk
// return NULL if run out of memory
T *get(int n) {
if (n > maxchunk) {
errorflag = 1;
return NULL;
}
ndatum += n;
nchunk++;
if (index+n <= pagesize) {
int start = index;
index += n;
return &page[start];
}
ipage++;
if (ipage == npage) {
allocate();
if (errorflag) return NULL;
}
page = pages[ipage];
index = n;
return &page[0];
}
// get ptr to location that can store maxchunk datums
// will return same ptr as previous call if vgot() not called
// return NULL if run out of memory
/** Get pointer to location that can store *maxchunk* items.
*
* This will return the same pointer as the previous call to
* this function unless vgot() is called afterwards to record
* how many items of the chunk were actually used.
*
* \return pointer to chunk of memory or null pointer if run out of memory */
T *vget() {
if (index+maxchunk <= pagesize) return &page[index];
@ -163,9 +60,14 @@ class MyPage {
return &page[index];
}
// increment by N = # of values stored in loc returned by vget()
// OK to not call if vget() ptr was not used
// error if N > maxchunk
/** Mark *N* items as used of the chunk reserved with a preceding call to vget().
*
* This will advance the internal pointer inside the current memory page.
* It is not necessary to call this function for *N* = 0, that is the reserved
* storage was not used. A following call to vget() will then reserve the
* same location again. It is an error if *N* > *maxchunk*.
*
* \param n Number of items used in previously reserved chunk */
void vgot(int n) {
if (n > maxchunk) errorflag = 1;
@ -174,25 +76,21 @@ class MyPage {
index += n;
}
// clear all pages, without freeing any memory
void reset();
void reset() {
ndatum = nchunk = 0;
index = ipage = 0;
page = pages[ipage];
/** Return total size of allocated pages
*
* \return total storage used in bytes */
double size() const {
return (double)npage*pagesize*sizeof(T);
}
// return total size of allocated pages
/** Return error status
*
* \return 0 if no error, 1 requested chunk size > maxchunk, 2 if malloc failed */
int size() const {
return npage*pagesize*sizeof(T);
}
// return error status
int status() const {
return errorflag;
}
int status() const { return errorflag; }
private:
T **pages; // list of allocated pages
@ -208,27 +106,8 @@ class MyPage {
int errorflag; // flag > 0 if error has occurred
// 1 = chunk size exceeded maxchunk
// 2 = memory allocation error
void allocate() {
npage += pagedelta;
pages = (T **) realloc(pages,npage*sizeof(T *));
if (!pages) {
errorflag = 2;
return;
}
for (int i = npage-pagedelta; i < npage; i++) {
#if defined(LAMMPS_MEMALIGN)
void *ptr;
if (posix_memalign(&ptr, LAMMPS_MEMALIGN, pagesize*sizeof(T)))
errorflag = 2;
pages[i] = (T *) ptr;
#else
pages[i] = (T *) malloc(pagesize*sizeof(T));
if (!pages[i]) errorflag = 2;
#endif
}
}
void allocate();
void deallocate();
};
}

228
src/my_pool_chunk.cpp Normal file
View File

@ -0,0 +1,228 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#include "my_pool_chunk.h"
#include <cstdlib>
#include <cstdio>
#if defined(LMP_USER_INTEL) && !defined(LAMMPS_MEMALIGN) && !defined(_WIN32)
#define LAMMPS_MEMALIGN 64
#endif
using namespace LAMMPS_NS;
/** \class LAMMPS_NS::MyPoolChunk
* \brief Templated class for storing chunks of datums in pages
*
* The size of the chunk may vary from call to call between the
* *minchunk* and *maxchunk* setting. Chunks may be returned
* to the pool for re-use. Chunks can be reserved in *nbin*
* different sizes between *minchunk* and *maxchunk*.
* The *chunksperpage* setting specifies how many chunks are stored
* on any page and the *pagedelta* setting determines how many
* pages are allocated in one go. Pages are never freed, so they
* can be re-used without re-allocation.
*
* \note
* This is a template class with explicit instantiation. If the class
* is used with a new data type a new explicit instantiation may need
* to be added at the end of the file ``src/my_pool_chunk.cpp`` to
* avoid symbol lookup errors. */
/** Create a class instance and set memory pool parameters
*
* \param user_minchunk Minimal chunk size
* \param user_maxchunk Maximal chunk size
* \param user_nbin Number of bins of different chunk sizes
* \param user_chunkperpage Number of chunks per page
* \param user_pagedelta Number of pages to allocate in one go */
template <class T>
MyPoolChunk<T>::MyPoolChunk(int user_minchunk, int user_maxchunk, int user_nbin,
int user_chunkperpage, int user_pagedelta) {
minchunk = user_minchunk;
maxchunk = user_maxchunk;
nbin = user_nbin;
chunkperpage = user_chunkperpage;
pagedelta = user_pagedelta;
errorflag = 0;
if (minchunk <= 0 || minchunk > maxchunk) errorflag = 1;
if (user_nbin <= 0 || chunkperpage <= 0 || pagedelta <= 0) errorflag = 1;
freehead = new int[nbin];
chunksize = new int[nbin];
if (!freehead || !chunksize) errorflag = 1;
if (errorflag) return;
// insure nbin*binsize spans minchunk to maxchunk inclusive
binsize = (maxchunk-minchunk+1) / nbin;
if (minchunk + nbin*binsize <= maxchunk) binsize++;
freelist = nullptr;
for (int ibin = 0; ibin < nbin; ibin++) {
freehead[ibin] = -1;
chunksize[ibin] = minchunk + (ibin+1)*binsize - 1;
if (chunksize[ibin] > maxchunk) chunksize[ibin] = maxchunk;
}
ndatum = nchunk = 0;
pages = nullptr;
whichbin = nullptr;
npage = 0;
}
/** Destroy class instance and free all allocated memory */
template <class T>
MyPoolChunk<T>::~MyPoolChunk() {
delete [] freehead;
delete [] chunksize;
if (npage) {
free(freelist);
for (int i = 0; i < npage; i++) free(pages[i]);
free(pages);
free(whichbin);
}
}
/** Return pointer/index of unused chunk of size maxchunk
*
* \param index Index of chunk in memory pool
* \return Pointer to requested chunk of storage */
template <class T>
T *MyPoolChunk<T>::get(int &index) {
int ibin = nbin-1;
if (freehead[ibin] < 0) {
allocate(ibin);
if (errorflag) {
index = -1;
return nullptr;
}
}
ndatum += maxchunk;
nchunk++;
index = freehead[ibin];
int ipage = index/chunkperpage;
int ientry = index % chunkperpage;
freehead[ibin] = freelist[index];
return &pages[ipage][ientry*chunksize[ibin]];
}
/** Return pointer/index of unused chunk of size N
*
* \param n Size of chunk
* \param index Index of chunk in memory pool
* \return Pointer to requested chunk of storage */
template <class T>
T *MyPoolChunk<T>::get(int n, int &index) {
if (n < minchunk || n > maxchunk) {
errorflag = 3;
index = -1;
return nullptr;
}
int ibin = (n-minchunk) / binsize;
if (freehead[ibin] < 0) {
allocate(ibin);
if (errorflag) {
index = -1;
return nullptr;
}
}
ndatum += n;
nchunk++;
index = freehead[ibin];
int ipage = index/chunkperpage;
int ientry = index % chunkperpage;
freehead[ibin] = freelist[index];
return &pages[ipage][ientry*chunksize[ibin]];
}
/** Put indexed chunk back into memory pool via free list
*
* \param index Memory chunk index returned by call to get() */
template <class T>
void MyPoolChunk<T>::put(int index) {
if (index < 0) return;
int ipage = index/chunkperpage;
int ibin = whichbin[ipage];
nchunk--;
ndatum -= chunksize[ibin];
freelist[index] = freehead[ibin];
freehead[ibin] = index;
}
template <class T>
void MyPoolChunk<T>::allocate(int ibin) {
int oldpage = npage;
npage += pagedelta;
freelist = (int *) realloc(freelist,npage*chunkperpage*sizeof(int));
pages = (T **) realloc(pages,npage*sizeof(T *));
whichbin = (int *) realloc(whichbin,npage*sizeof(int));
if (!freelist || !pages) {
errorflag = 2;
return;
}
// allocate pages with appropriate chunksize for ibin
for (int i = oldpage; i < npage; i++) {
whichbin[i] = ibin;
#if defined(LAMMPS_MEMALIGN)
void *ptr;
if (posix_memalign(&ptr, LAMMPS_MEMALIGN,
chunkperpage*chunksize[ibin]*sizeof(T)))
errorflag = 2;
pages[i] = (T *) ptr;
#else
pages[i] = (T *) malloc(chunkperpage*chunksize[ibin]*sizeof(T));
if (!pages[i]) errorflag = 2;
#endif
}
// reset free list for unused chunks on new pages
freehead[ibin] = oldpage*chunkperpage;
for (int i = freehead[ibin]; i < npage*chunkperpage; i++) freelist[i] = i+1;
freelist[npage*chunkperpage-1] = -1;
}
/** Return total size of allocated pages
*
* \return total storage used in bytes */
template <class T>
double MyPoolChunk<T>::size() const {
double bytes = npage*chunkperpage*sizeof(int);
bytes += npage*sizeof(T *);
bytes += npage*sizeof(int);
for (int i=0; i < npage; ++i)
bytes += chunkperpage*chunksize[i]*sizeof(T);
return bytes;
}
// explicit instantiations
namespace LAMMPS_NS {
template class MyPoolChunk<int>;
template class MyPoolChunk<double>;
}

View File

@ -9,46 +9,11 @@
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
MyPoolChunk = templated class for storing chunks of datums in pages
chunks can be returned to pool for reuse
chunks come in nbin different fixed sizes so can reuse
replaces many small mallocs with a few large mallocs
pages are never freed, so can reuse w/out reallocs
usage:
continuously get() and put() chunks as needed
NOTE: could add a clear() if retain info on mapping of pages to bins
inputs:
template T = one datum, e.g. int, double, struct
minchunk = min # of datums in one chunk, def = 1
maxchunk = max # of datums in one chunk, def = 1
nbin = # of bins between minchunk and maxchunk
chunkperpage = # of chunks in one page, def = 1024
pagedelta = # of pages to allocate at a time, def = 1
methods:
T *get(index) = return ptr/index to unused chunk of size maxchunk
T *get(N,index) = return ptr/index to unused chunk of size N
minchunk <= N <= maxchunk required
put(index) = return indexed chunk to pool (same index returned by get)
int size() = return total size of allocated pages in bytes
public variables:
ndatum = total # of stored datums
nchunk = total # of stored chunks
size = total size of all allocated pages in daums
errorflag = flag for various error conditions
------------------------------------------------------------------------- */
------------------------------------------------------------------------- */
#ifndef LAMMPS_MY_POOL_CHUNK_H
#define LAMMPS_MY_POOL_CHUNK_H
#if defined(LMP_USER_INTEL) && !defined(LAMMPS_MEMALIGN) && !defined(_WIN32)
#define LAMMPS_MEMALIGN 64
#endif
#include <cstdlib>
namespace LAMMPS_NS {
template<class T>
@ -56,113 +21,36 @@ class MyPoolChunk {
public:
int ndatum; // total # of stored datums
int nchunk; // total # of stored chunks
int size; // total size of all allocated pages in datums
int errorflag; // flag > 1 if error has occurred
// 1 = invalid inputs
// 2 = memory allocation error
// 3 = chunk size exceeded maxchunk
MyPoolChunk(int user_minchunk = 1, int user_maxchunk = 1, int user_nbin = 1,
int user_chunkperpage = 1024, int user_pagedelta = 1) {
minchunk = user_minchunk;
maxchunk = user_maxchunk;
nbin = user_nbin;
chunkperpage = user_chunkperpage;
pagedelta = user_pagedelta;
errorflag = 0;
if (minchunk <= 0 || minchunk > maxchunk) errorflag = 1;
if (user_nbin <= 0 || chunkperpage <= 0 || pagedelta <= 0) errorflag = 1;
freehead = new int[nbin];
chunksize = new int[nbin];
if (!freehead || !chunksize) errorflag = 1;
if (errorflag) return;
// insure nbin*binsize spans minchunk to maxchunk inclusive
binsize = (maxchunk-minchunk+1) / nbin;
if (minchunk + nbin*binsize <= maxchunk) binsize++;
freelist = NULL;
for (int ibin = 0; ibin < nbin; ibin++) {
freehead[ibin] = -1;
chunksize[ibin] = minchunk + (ibin+1)*binsize - 1;
if (chunksize[ibin] > maxchunk) chunksize[ibin] = maxchunk;
}
ndatum = nchunk = size = 0;
pages = NULL;
whichbin = NULL;
npage = 0;
}
int user_chunkperpage = 1024, int user_pagedelta = 1);
// free all allocated memory
~MyPoolChunk() {
delete [] freehead;
delete [] chunksize;
if (npage) {
free(freelist);
for (int i = 0; i < npage; i++) free(pages[i]);
free(pages);
free(whichbin);
}
}
~MyPoolChunk();
// return pointer/index of unused chunk of size maxchunk
T *get(int &index) {
int ibin = nbin-1;
if (freehead[ibin] < 0) {
allocate(ibin);
if (errorflag) return NULL;
}
ndatum += maxchunk;
nchunk++;
index = freehead[ibin];
int ipage = index/chunkperpage;
int ientry = index % chunkperpage;
freehead[ibin] = freelist[index];
return &pages[ipage][ientry*chunksize[ibin]];
}
T *get(int &index);
// return pointer/index of unused chunk of size N
T *get(int n, int &index) {
if (n < minchunk || n > maxchunk) {
errorflag = 3;
return NULL;
}
int ibin = (n-minchunk) / binsize;
if (freehead[ibin] < 0) {
allocate(ibin);
if (errorflag) return NULL;
}
ndatum += n;
nchunk++;
index = freehead[ibin];
int ipage = index/chunkperpage;
int ientry = index % chunkperpage;
freehead[ibin] = freelist[index];
return &pages[ipage][ientry*chunksize[ibin]];
}
T *get(int n, int &index);
// return indexed chunk to pool via free list
// index = -1 if no allocated chunk
void put(int index) {
if (index < 0) return;
int ipage = index/chunkperpage;
int ibin = whichbin[ipage];
nchunk--;
ndatum -= chunksize[ibin];
freelist[index] = freehead[ibin];
freehead[ibin] = index;
}
void put(int index);
// total memory used in bytes
double size() const;
/** Return error status
*
* \return 0 if no error, 1 if invalid input, 2 if malloc() failed, 3 if chunk > maxchunk */
int status() const { return errorflag; }
private:
int minchunk; // min # of datums per chunk
@ -171,6 +59,10 @@ class MyPoolChunk {
int chunkperpage; // # of chunks on every page, regardless of which bin
int pagedelta; // # of pages to allocate at once, default = 1
int binsize; // delta in chunk sizes between adjacent bins
int errorflag; // flag > 0 if error has occurred
// 1 = invalid inputs
// 2 = memory allocation error
// 3 = chunk size exceeded maxchunk
T **pages; // list of allocated pages
int *whichbin; // which bin each page belongs to
@ -179,42 +71,7 @@ class MyPoolChunk {
int *freehead; // index of first unused chunk in each bin
int *chunksize; // size of chunks in each bin
void allocate(int ibin) {
int oldpage = npage;
npage += pagedelta;
freelist = (int *) realloc(freelist,npage*chunkperpage*sizeof(int));
pages = (T **) realloc(pages,npage*sizeof(T *));
whichbin = (int *) realloc(whichbin,npage*sizeof(int));
if (!freelist || !pages) {
errorflag = 2;
return;
}
// allocate pages with appropriate chunksize for ibin
for (int i = oldpage; i < npage; i++) {
whichbin[i] = ibin;
#if defined(LAMMPS_MEMALIGN)
void *ptr;
if (posix_memalign(&ptr, LAMMPS_MEMALIGN,
chunkperpage*chunksize[ibin]*sizeof(T)))
errorflag = 2;
pages[i] = (T *) ptr;
#else
pages[i] = (T *) malloc(chunkperpage*chunksize[ibin]*sizeof(T));
size += chunkperpage*chunksize[ibin];
if (!pages[i]) errorflag = 2;
#endif
}
// reset free list for unused chunks on new pages
freehead[ibin] = oldpage*chunkperpage;
for (int i = freehead[ibin]; i < npage*chunkperpage; i++) freelist[i] = i+1;
freelist[npage*chunkperpage-1] = -1;
}
void allocate(int ibin);
};
}
#endif

View File

@ -34,7 +34,6 @@ using namespace LAMMPS_NS;
ResetIDs::AtomRvous *ResetIDs::sortrvous;
static int compare_coords(const void *, const void *);
#else
#include "mergesort.h"
// prototype for non-class function
static int compare_coords(const int, const int, void *);
#endif
@ -509,7 +508,7 @@ int ResetIDs::sort_bins(int n, char *inbuf,
sortrvous = in;
qsort(order,count[ibin],sizeof(int),compare_coords);
#else
merge_sort(order,count[ibin],(void *) in,compare_coords);
utils::merge_sort(order,count[ibin],(void *) in,compare_coords);
#endif
head[ibin] = last[ibin] = -1;

View File

@ -71,6 +71,16 @@ extern "C"
static int re_match(const char *text, const char *pattern);
}
////////////////////////////////////////////////////////////////////////
// Merge sort support functions
static void do_merge(int *idx, int *buf, int llo, int lhi, int rlo, int rhi,
void *ptr, int (*comp)(int, int, void *));
static void insertion_sort(int *index, int num, void *ptr,
int (*comp)(int, int, void*));
////////////////////////////////////////////////////////////////////////
using namespace LAMMPS_NS;
/** More flexible and specific matching of a string against a pattern.
@ -1011,6 +1021,113 @@ int utils::date2num(const std::string &date)
return num;
}
/* ----------------------------------------------------------------------
* Merge sort part 1: Loop over sublists doubling in size with each iteration.
* Pre-sort small sublists with insertion sort for better overall performance.
------------------------------------------------------------------------- */
void utils::merge_sort(int *index, int num, void *ptr,
int (*comp)(int, int, void *))
{
if (num < 2) return;
int chunk,i,j;
// do insertion sort on chunks of up to 64 elements
chunk = 64;
for (i=0; i < num; i += chunk) {
j = (i+chunk > num) ? num-i : chunk;
insertion_sort(index+i,j,ptr,comp);
}
// already done?
if (chunk >= num) return;
// continue with merge sort on the pre-sorted chunks.
// we need an extra buffer for temporary storage and two
// pointers to operate on, so we can swap the pointers
// rather than copying to the hold buffer in each pass
int *buf = new int[num];
int *dest = index;
int *hold = buf;
while (chunk < num) {
int m;
// swap hold and destination buffer
int *tmp = dest; dest = hold; hold = tmp;
// merge from hold array to destination array
for (i=0; i < num-1; i += 2*chunk) {
j = i + 2*chunk;
if (j > num) j=num;
m = i+chunk;
if (m > num) m=num;
do_merge(dest,hold,i,m,m,j,ptr,comp);
}
// copy all indices not handled by the chunked merge sort loop
for ( ; i < num ; i++ ) dest[i] = hold[i];
chunk *= 2;
}
// if the final sorted data is in buf, copy back to index
if (dest == buf) memcpy(index,buf,sizeof(int)*num);
delete[] buf;
}
/* ------------------------------------------------------------------ */
/* ----------------------------------------------------------------------
* Merge sort part 2: Insertion sort for pre-sorting of small chunks
------------------------------------------------------------------------- */
void insertion_sort(int *index, int num, void *ptr,
int (*comp)(int, int, void*))
{
if (num < 2) return;
for (int i=1; i < num; ++i) {
int tmp = index[i];
for (int j=i-1; j >= 0; --j) {
if ((*comp)(index[j],tmp,ptr) > 0) {
index[j+1] = index[j];
} else {
index[j+1] = tmp;
break;
}
if (j == 0) index[0] = tmp;
}
}
}
/* ----------------------------------------------------------------------
* Merge sort part 3: Merge two sublists
------------------------------------------------------------------------- */
static void do_merge(int *idx, int *buf, int llo, int lhi, int rlo, int rhi,
void *ptr, int (*comp)(int, int, void *))
{
int i = llo;
int l = llo;
int r = rlo;
while ((l < lhi) && (r < rhi)) {
if ((*comp)(buf[l],buf[r],ptr) < 0)
idx[i++] = buf[l++];
else idx[i++] = buf[r++];
}
while (l < lhi) idx[i++] = buf[l++];
while (r < rhi) idx[i++] = buf[r++];
}
/* ------------------------------------------------------------------ */
extern "C" {

View File

@ -33,23 +33,23 @@ namespace LAMMPS_NS {
*
* \param text the text to be matched against the pattern
* \param pattern the search pattern, which may contain regexp markers
* \return true if the pattern matches, false if not
*/
* \return true if the pattern matches, false if not */
bool strmatch(const std::string &text, const std::string &pattern);
/** Send message to screen and logfile, if available
*
* \param lmp pointer to LAMMPS class instance
* \param mesg message to be printed
*/
* \param mesg message to be printed */
void logmesg(LAMMPS *lmp, const std::string &mesg);
/** return a string representing the current system error status
*
* This is a wrapper around calling strerror(errno).
*
* \return error string
*/
* \return error string */
std::string getsyserror();
/** safe wrapper around fgets() which aborts on errors
@ -61,8 +61,8 @@ namespace LAMMPS_NS {
* \param size size of buffer s (max number of bytes read by fgets())
* \param fp file pointer used by fgets()
* \param filename file name associated with fp (may be NULL; then LAMMPS will try to detect)
* \param error pointer to Error class instance (for abort)
*/
* \param error pointer to Error class instance (for abort) */
void sfgets(const char *srcname, int srcline, char *s, int size,
FILE *fp, const char *filename, Error *error);
@ -76,8 +76,8 @@ namespace LAMMPS_NS {
* \param num number of data elements read by fread()
* \param fp file pointer used by fread()
* \param filename file name associated with fp (may be NULL; then LAMMPS will try to detect)
* \param error pointer to Error class instance (for abort)
*/
* \param error pointer to Error class instance (for abort) */
void sfread(const char *srcname, int srcline, void *s, size_t size,
size_t num, FILE *fp, const char *filename, Error *error);
@ -86,8 +86,8 @@ namespace LAMMPS_NS {
* \param style type of style that is to be checked for
* \param name name of style that was not found
* \param lmp pointer to top-level LAMMPS class instance
* \return string usable for error messages
*/
* \return string usable for error messages */
std::string check_packages_for_style(const std::string &style,
const std::string &name, LAMMPS *lmp);
@ -112,8 +112,8 @@ namespace LAMMPS_NS {
* \param str string to be converted to number
* \param do_abort determines whether to call Error::one() or Error::all()
* \param lmp pointer to top-level LAMMPS class instance
* \return integer number (regular int)
*/
* \return integer number (regular int) */
int inumeric(const char *file, int line, const char *str,
bool do_abort, LAMMPS *lmp);
@ -125,8 +125,8 @@ namespace LAMMPS_NS {
* \param str string to be converted to number
* \param do_abort determines whether to call Error::one() or Error::all()
* \param lmp pointer to top-level LAMMPS class instance
* \return integer number (bigint)
*/
* \return integer number (bigint) */
bigint bnumeric(const char *file, int line, const char *str,
bool do_abort, LAMMPS *lmp);
@ -162,6 +162,7 @@ namespace LAMMPS_NS {
* \param nlo lower bound
* \param nhi upper bound
* \param error pointer to Error class for out-of-bounds messages */
template <typename TYPE>
void bounds(const char *file, int line, const std::string &str,
bigint nmin, bigint nmax, TYPE &nlo, TYPE &nhi, Error *error);
@ -197,45 +198,45 @@ namespace LAMMPS_NS {
/** Trim leading and trailing whitespace. Like TRIM() in Fortran.
*
* \param line string that should be trimmed
* \return new string without whitespace (string)
*/
* \return new string without whitespace (string) */
std::string trim(const std::string &line);
/** Return string with anything from '#' onward removed
*
* \param line string that should be trimmed
* \return new string without comment (string)
*/
* \return new string without comment (string) */
std::string trim_comment(const std::string &line);
/** Count words in string with custom choice of separating characters
*
* \param text string that should be searched
* \param separators string containing characters that will be treated as whitespace
* \return number of words found
*/
* \return number of words found */
size_t count_words(const std::string &text, const std::string &separators);
/** Count words in string, ignore any whitespace matching " \t\r\n\f"
*
* \param text string that should be searched
* \return number of words found
*/
* \return number of words found */
size_t count_words(const std::string &text);
/** Count words in C-string, ignore any whitespace matching " \t\r\n\f"
*
* \param text string that should be searched
* \return number of words found
*/
* \return number of words found */
size_t count_words(const char *text);
/** Count words in a single line, trim anything from '#' onward
*
* \param text string that should be trimmed and searched
* \param separators string containing characters that will be treated as whitespace
* \return number of words found
*/
* \return number of words found */
size_t trim_and_count_words(const std::string &text, const std::string &separators = " \t\r\n\f");
/** Take text and split into non-whitespace words.
@ -247,22 +248,22 @@ namespace LAMMPS_NS {
* Use a tokenizer class for that.
*
* \param text string that should be split
* \return STL vector with the words
*/
* \return STL vector with the words */
std::vector<std::string> split_words(const std::string &text);
/** Check if string can be converted to valid integer
*
* \param str string that should be checked
* \return true, if string contains valid integer, false otherwise
*/
* \return true, if string contains valid a integer, false otherwise */
bool is_integer(const std::string &str);
/** Check if string can be converted to valid floating-point number
*
* \param str string that should be checked
* \return true, if string contains valid floating-point number, false otherwise
*/
* \return true, if string contains valid number, false otherwise */
bool is_double(const std::string &str);
/** Try to detect pathname from FILE pointer.
@ -272,55 +273,60 @@ namespace LAMMPS_NS {
* \param buf storage buffer for pathname. output will be truncated if not large enough
* \param len size of storage buffer. output will be truncated to this length - 1
* \param fp FILE pointer struct from STDIO library for which we want to detect the name
* \return pointer to the storage buffer, i.e. buf
*/
* \return pointer to the storage buffer, i.e. buf */
const char *guesspath(char *buf, int len, FILE *fp);
/** Strip off leading part of path, return just the filename
*
* \param path file path
* \return file name
*/
* \return file name */
std::string path_basename(const std::string &path);
/**
* \brief Join two paths
/** Join two pathname segments
*
* This uses the forward slash '/' character unless LAMMPS is compiled
* for Windows where it used the equivalent backward slash '\\'.
*
* \param a first path
* \param b second path
* \return combined path
*/
* \return combined path */
std::string path_join(const std::string &a, const std::string &b);
/**
* \brief Check if file exists and is readable
/** Check if file exists and is readable
*
* \param path file path
* \return true if file exists and is readable
*/
* \return true if file exists and is readable */
bool file_is_readable(const std::string &path);
/** Determine full path of potential file. If file is not found in current directory,
* search directories listed in LAMMPS_POTENTIALS environment variable
*
* \param path file path
* \return full path to potential file
*/
* \return full path to potential file */
std::string get_potential_file_path(const std::string &path);
/** Read potential file and return DATE field if it is present
*
* \param path file path
* \param potential_name name of potential that is being read
* \return DATE field if present
*/
std::string get_potential_date(const std::string &path, const std::string &potential_name);
* \return DATE field if present */
std::string get_potential_date(const std::string &path,
const std::string &potential_name);
/** Read potential file and return UNITS field if it is present
*
* \param path file path
* \param potential_name name of potential that is being read
* \return UNITS field if present
*/
std::string get_potential_units(const std::string &path, const std::string &potential_name);
* \return UNITS field if present */
std::string get_potential_units(const std::string &path,
const std::string &potential_name);
enum { NOCONVERT = 0, METAL2REAL = 1, REAL2METAL = 1<<1 };
enum { UNKNOWN = 0, ENERGY };
@ -328,16 +334,15 @@ namespace LAMMPS_NS {
/** Return bitmask of available conversion factors for a given property
*
* \param property property to be converted
* \return bitmask indicating available conversions
*/
* \return bitmask indicating available conversions */
int get_supported_conversions(const int property);
/** Return unit conversion factor for given property and selected from/to units
*
* \param property property to be converted
* \param conversion constant indicating the conversion
* \return conversion factor
*/
* \return conversion factor */
double get_conversion_factor(const int property, const int conversion);
/** Open a potential file as specified by *name*
@ -368,8 +373,8 @@ namespace LAMMPS_NS {
* The strings "off" and "unlimited" result in -1
*
* \param timespec a string in the following format: ([[HH:]MM:]SS)
* \return total in seconds
*/
* \return total in seconds */
double timespec2seconds(const std::string &timespec);
/** Convert a LAMMPS version date to a number
@ -386,9 +391,26 @@ namespace LAMMPS_NS {
* No check is made whether the date is valid.
*
* \param date string in the format (Day Month Year)
* \return date code
*/
* \return date code */
int date2num(const std::string &date);
/** Custom merge sort implementation
*
* This function provides a custom upward hybrid merge sort
* implementation with support to pass an opaque pointer to
* the comparison function, e.g. for access to class members.
* This avoids having to use global variables. For improved
* performance, it uses an in-place insertion sort on initial
* chunks of up to 64 elements and switches to merge sort from
* then on.
*
* \param index Array with indices to be sorted
* \param num Length of the index array
* \param ptr Pointer to opaque object passed to comparison function
* \param comp Pointer to comparison function */
void merge_sort(int *index, int num, void *ptr,
int (*comp)(int, int, void *));
}
}

View File

@ -164,7 +164,7 @@ TEST(lammps_open_fortran, no_args) {
MPI_Comm_split(MPI_COMM_WORLD, 0, 1, &mycomm);
int fcomm = MPI_Comm_c2f(mycomm);
::testing::internal::CaptureStdout();
void *handle = lammps_open_fortran(0, NULL, fcomm, NULL);
void *handle = lammps_open_fortran(0, NULL, fcomm);
std::string output = ::testing::internal::GetCapturedStdout();
EXPECT_STREQ(output.substr(0,6).c_str(),"LAMMPS");
LAMMPS_NS::LAMMPS *lmp = (LAMMPS_NS::LAMMPS *)handle;

View File

@ -2,6 +2,10 @@ add_executable(test_tokenizer test_tokenizer.cpp)
target_link_libraries(test_tokenizer PRIVATE lammps GTest::GMockMain GTest::GMock GTest::GTest)
add_test(Tokenizer test_tokenizer)
add_executable(test_mempool test_mempool.cpp)
target_link_libraries(test_mempool PRIVATE lammps GTest::GMockMain GTest::GMock GTest::GTest)
add_test(MemPool test_mempool)
add_executable(test_utils test_utils.cpp)
target_link_libraries(test_utils PRIVATE lammps GTest::GMockMain GTest::GMock GTest::GTest)
add_test(Utils test_utils)

View File

@ -0,0 +1,347 @@
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#include "lmptype.h"
#include "my_page.h"
#include "my_pool_chunk.h"
#include "gmock/gmock.h"
#include "gtest/gtest.h"
using namespace LAMMPS_NS;
TEST(MyPage, int) {
MyPage<int> p;
// default init. maxchunk=1, pagesize=1024
int rv = p.init();
ASSERT_EQ(rv,0);
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
int *iptr = p.vget();
// second call to vget() should give same pointer without vgot()
ASSERT_EQ(iptr,p.vget());
p.vgot(1);
++iptr;
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,1);
ASSERT_EQ(p.nchunk,1);
ASSERT_EQ(iptr,p.vget());
// use too large chunk size
p.vgot(2);
ASSERT_EQ(1,p.status());
p.reset();
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
p.vgot(1);
++iptr;
ASSERT_EQ(iptr,p.get());
++iptr;
ASSERT_EQ(iptr,p.get(1));
ASSERT_EQ(p.ndatum,3);
ASSERT_EQ(p.nchunk,3);
// restart with custom init. maxchunk=16, pagesize=256
rv = p.init(16,64,2);
ASSERT_EQ(rv,0);
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
// second call to vget() should give same pointer without vgot()
ASSERT_EQ(iptr,p.vget());
p.vgot(16);
iptr += 16;
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,16);
ASSERT_EQ(p.nchunk,1);
// use too large chunk size
ASSERT_EQ(iptr,p.vget());
p.vgot(32);
ASSERT_EQ(1,p.status());
p.reset();
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
p.vgot(16);
iptr = p.vget();
p.vgot(4);
iptr += 4;
ASSERT_EQ(iptr,p.get());
++iptr;
ASSERT_EQ(iptr,p.get(16));
ASSERT_DOUBLE_EQ(p.size(),(double)sizeof(int)*128.0);
ASSERT_EQ(p.ndatum,37);
ASSERT_EQ(p.nchunk,4);
p.get(16);
p.get(16);
// allocation on the same page
iptr = p.get(16);
iptr += 16;
ASSERT_EQ(iptr,p.get(16));
// allocation on different pages
p.get(16);
iptr += 16;
ASSERT_NE(iptr,p.get(16));
ASSERT_DOUBLE_EQ(p.size(),(double)sizeof(int)*256.0);
ASSERT_EQ(p.ndatum,133);
ASSERT_EQ(p.nchunk,10);
}
TEST(MyPage, double) {
MyPage<double> p;
// default init. maxchunk=1, pagesize=1024
int rv = p.init();
ASSERT_EQ(rv,0);
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
double *iptr = p.vget();
// second call to vget() should give same pointer without vgot()
ASSERT_EQ(iptr,p.vget());
p.vgot(1);
++iptr;
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,1);
ASSERT_EQ(p.nchunk,1);
ASSERT_EQ(iptr,p.vget());
// use too large chunk size
p.vgot(2);
ASSERT_EQ(1,p.status());
p.reset();
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
p.vgot(1);
++iptr;
ASSERT_EQ(iptr,p.get());
++iptr;
ASSERT_EQ(iptr,p.get(1));
ASSERT_EQ(p.ndatum,3);
ASSERT_EQ(p.nchunk,3);
// restart with custom init. maxchunk=16, pagesize=256
rv = p.init(16,64,2);
ASSERT_EQ(rv,0);
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
// second call to vget() should give same pointer without vgot()
ASSERT_EQ(iptr,p.vget());
p.vgot(16);
iptr += 16;
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,16);
ASSERT_EQ(p.nchunk,1);
// use too large chunk size
ASSERT_EQ(iptr,p.vget());
p.vgot(32);
ASSERT_EQ(1,p.status());
p.reset();
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
p.vgot(16);
iptr = p.vget();
p.vgot(4);
iptr += 4;
ASSERT_EQ(iptr,p.get());
++iptr;
ASSERT_EQ(iptr,p.get(16));
ASSERT_DOUBLE_EQ(p.size(),(double)sizeof(double)*128.0);
ASSERT_EQ(p.ndatum,37);
ASSERT_EQ(p.nchunk,4);
p.get(16);
p.get(16);
// allocation on the same page
iptr = p.get(16);
iptr += 16;
ASSERT_EQ(iptr,p.get(16));
// allocation on different pages
p.get(16);
iptr += 16;
ASSERT_NE(iptr,p.get(16));
ASSERT_DOUBLE_EQ(p.size(),(double)sizeof(double)*256.0);
ASSERT_EQ(p.ndatum,133);
ASSERT_EQ(p.nchunk,10);
}
TEST(MyPage, bigint) {
MyPage<bigint> p;
// default init. maxchunk=1, pagesize=1024
int rv = p.init();
ASSERT_EQ(rv,0);
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
bigint *iptr = p.vget();
// second call to vget() should give same pointer without vgot()
ASSERT_EQ(iptr,p.vget());
p.vgot(1);
++iptr;
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,1);
ASSERT_EQ(p.nchunk,1);
ASSERT_EQ(iptr,p.vget());
// use too large chunk size
p.vgot(2);
ASSERT_EQ(1,p.status());
p.reset();
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
p.vgot(1);
++iptr;
ASSERT_EQ(iptr,p.get());
++iptr;
ASSERT_EQ(iptr,p.get(1));
ASSERT_EQ(p.ndatum,3);
ASSERT_EQ(p.nchunk,3);
// restart with custom init. maxchunk=16, pagesize=256
rv = p.init(16,64,2);
ASSERT_EQ(rv,0);
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
// second call to vget() should give same pointer without vgot()
ASSERT_EQ(iptr,p.vget());
p.vgot(16);
iptr += 16;
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,16);
ASSERT_EQ(p.nchunk,1);
// use too large chunk size
ASSERT_EQ(iptr,p.vget());
p.vgot(32);
ASSERT_EQ(1,p.status());
p.reset();
ASSERT_EQ(0,p.status());
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
iptr = p.vget();
p.vgot(16);
iptr = p.vget();
p.vgot(4);
iptr += 4;
ASSERT_EQ(iptr,p.get());
++iptr;
ASSERT_EQ(iptr,p.get(16));
ASSERT_DOUBLE_EQ(p.size(),(double)sizeof(bigint)*128.0);
ASSERT_EQ(p.ndatum,37);
ASSERT_EQ(p.nchunk,4);
p.get(16);
p.get(16);
// allocation on the same page
iptr = p.get(16);
iptr += 16;
ASSERT_EQ(iptr,p.get(16));
// allocation on different pages
p.get(16);
iptr += 16;
ASSERT_NE(iptr,p.get(16));
ASSERT_DOUBLE_EQ(p.size(),(double)sizeof(bigint)*256.0);
ASSERT_EQ(p.ndatum,133);
ASSERT_EQ(p.nchunk,10);
}
TEST(MyPoolChunk, int) {
// defaults to minchunk=1, maxchunk=1, nbin=1, chunksperpage=1024, pagedelta=1
MyPoolChunk<int> p;
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
ASSERT_EQ(p.size(),0.0);
int idx=~0x0000;
int *iptr = p.get(idx);
ASSERT_NE(iptr,nullptr);
ASSERT_EQ(idx,0);
iptr = p.get(1,idx);
ASSERT_NE(iptr,nullptr);
ASSERT_EQ(idx,1);
// we have only one page allocated
ASSERT_EQ(p.size(),1024*sizeof(int)+1024*sizeof(int)+sizeof(void *)+sizeof(int));
ASSERT_EQ(p.ndatum,2);
ASSERT_EQ(p.nchunk,2);
p.put(0);
ASSERT_EQ(p.ndatum,1);
ASSERT_EQ(p.nchunk,1);
iptr = p.get(2,idx);
ASSERT_EQ(iptr,nullptr);
ASSERT_EQ(p.status(),3);
ASSERT_EQ(p.ndatum,1);
ASSERT_EQ(p.nchunk,1);
}
TEST(MyPoolChunk, double) {
// defaults to minchunk=1, maxchunk=1, nbin=1, chunksperpage=1024, pagedelta=1
MyPoolChunk<double> p;
ASSERT_EQ(p.ndatum,0);
ASSERT_EQ(p.nchunk,0);
ASSERT_EQ(p.size(),0.0);
int idx=~0x0000;
double *dptr = p.get(idx);
ASSERT_NE(dptr,nullptr);
ASSERT_EQ(idx,0);
dptr = p.get(1,idx);
ASSERT_NE(dptr,nullptr);
ASSERT_EQ(idx,1);
// we have only one page allocated
ASSERT_EQ(p.size(),1024*sizeof(int)+1024*sizeof(double)+sizeof(void *)+sizeof(int));
p.put(0);
ASSERT_EQ(p.ndatum,1);
ASSERT_EQ(p.nchunk,1);
dptr = p.get(2,idx);
ASSERT_EQ(dptr,nullptr);
ASSERT_EQ(p.status(),3);
ASSERT_EQ(p.ndatum,1);
ASSERT_EQ(p.nchunk,1);
}