first try at implementing lammps_extract_atom_size()

doc/src/Library_atoms.rst

@@ -4,6 +4,7 @@ Per-atom properties
 This section documents the following functions:
 
 - :cpp:func:`lammps_extract_atom_datatype`
+- :cpp:func:`lammps_extract_atom_size`
 - :cpp:func:`lammps_extract_atom`
 
 -----------------------
@@ -13,6 +14,11 @@ This section documents the following functions:
 
 -----------------------
 
+.. doxygenfunction:: lammps_extract_atom_size
+   :project: progguide
+
+-----------------------
+
 .. doxygenfunction:: lammps_extract_atom
    :project: progguide
 

examples/COUPLE/plugin/liblammpsplugin.c

@@ -105,6 +105,7 @@ liblammpsplugin_t *liblammpsplugin_load(const char *lib)
   ADDSYM(map_atom);
 
   ADDSYM(extract_atom_datatype);
+  ADDSYM(extract_atom_size);
   ADDSYM(extract_atom);
 
   ADDSYM(extract_compute);

examples/COUPLE/plugin/liblammpsplugin.h

@@ -151,6 +151,7 @@ struct _liblammpsplugin {
   int (*map_atom)(void *, const void *);
 
   int (*extract_atom_datatype)(void *, const char *);
+  int (*extract_atom_size)(void *, const char *, int);
   void *(*extract_atom)(void *, const char *);
 
   void *(*extract_compute)(void *, const char *, int, int);

fortran/lammps.f90

@@ -542,6 +542,14 @@ MODULE LIBLAMMPS
       INTEGER(c_int) :: lammps_extract_atom_datatype
     END FUNCTION lammps_extract_atom_datatype
 
+    FUNCTION lammps_extract_atom_size(handle, name, dtype) BIND(C)
+      IMPORT :: c_ptr, c_int
+      IMPLICIT NONE
+      TYPE(c_ptr), INTENT(IN), VALUE :: handle, name
+      INTEGER(c_int), INTENT(IN), VALUE :: dtype
+      INTEGER(c_int) :: lammps_extract_atom_size
+    END FUNCTION lammps_extract_atom_size
+
     FUNCTION lammps_extract_atom(handle, name) BIND(C)
       IMPORT :: c_ptr
       IMPLICIT NONE
@@ -1461,21 +1469,12 @@ CONTAINS
     ntypes = lmp_extract_setting(self, 'ntypes')
     Cname = f2c_string(name)
     datatype = lammps_extract_atom_datatype(self%handle, Cname)
+    ! Fortran and C/C++ have rows and columns switched
+    ncols = lammps_extract_atom_size(self%handle, Cname, LMP_SIZE_ROWS)
+    nrows = lammps_extract_atom_size(self%handle, Cname, LMP_SIZE_COLS)
     Cptr = lammps_extract_atom(self%handle, Cname)
     CALL lammps_free(Cname)
 
-    SELECT CASE (name)
-      CASE ('mass')
-        ncols = ntypes + 1
-        nrows = 1
-      CASE ('x','v','f','mu','omega','torque','angmom')
-        ncols = nmax
-        nrows = 3
-      CASE DEFAULT
-        ncols = nmax
-        nrows = 1
-    END SELECT
-
     peratom_data%lammps_instance => self
     SELECT CASE (datatype)
       CASE (LAMMPS_INT)
@@ -1486,6 +1485,7 @@ CONTAINS
         CALL C_F_POINTER(Cptr, peratom_data%i64_vec, [ncols])
       CASE (LAMMPS_DOUBLE)
         peratom_data%datatype = DATA_DOUBLE_1D
+        ! The mass array is allocated from 0, but only used from 1. We also want to use it from 1.
         IF (name == 'mass') THEN
           CALL C_F_POINTER(Cptr, dummy, [ncols])
           peratom_data%r64_vec(0:) => dummy

python/lammps/core.py

@@ -318,6 +318,8 @@ class lammps(object):
     self.lib.lammps_extract_atom.argtypes = [c_void_p, c_char_p]
     self.lib.lammps_extract_atom_datatype.argtypes = [c_void_p, c_char_p]
     self.lib.lammps_extract_atom_datatype.restype = c_int
+    self.lib.lammps_extract_atom_size.argtypes = [c_void_p, c_char_p, c_int]
+    self.lib.lammps_extract_atom_size.restype = c_int
 
     self.lib.lammps_extract_fix.argtypes = [c_void_p, c_char_p, c_int, c_int, c_int, c_int]
 
@@ -1070,31 +1072,59 @@ class lammps(object):
     else: return None
     return self.lib.lammps_extract_atom_datatype(self.lmp, newname)
 
+  # -------------------------------------------------------------------------
+  # extract per-atom info datatype
+
+  def extract_atom_size(self, name, dtype):
+    """Retrieve per-atom property dimensions from LAMMPS
+
+    This is a wrapper around the :cpp:func:`lammps_extract_atom_size`
+    function of the C-library interface. Its documentation includes a
+    list of the supported keywords.
+    This function returns ``None`` if the keyword is not
+    recognized. Otherwise it will return an integer value with the size
+    of the per-atom vector or array.  If *name* corresponds to a per-atom
+    array, the *dtype* keyword must be either LMP_SIZE_ROWS or LMP_SIZE_COLS
+    from the :ref:`type <py_type_constants>` constants defined in the
+    :py:mod:`lammps` module.  The return value is the requested size.
+    If *name* corresponds to a per-atom vector the *dtype* keyword is ignored.
+
+    :param name: name of the property
+    :type name:  string
+    :param type: either LMP_SIZE_ROWS or LMP_SIZE_COLS for arrays, otherwise ignored
+    :type type:  int
+    :return: data type of per-atom property (see :ref:`py_datatype_constants`)
+    :rtype: int
+    """
+    if name: newname = name.encode()
+    else: return None
+    return self.lib.lammps_extract_atom_size(self.lmp, newname, dtype)
+
   # -------------------------------------------------------------------------
   # extract per-atom info
 
   def extract_atom(self, name, dtype=LAMMPS_AUTODETECT):
     """Retrieve per-atom properties from LAMMPS
 
-    This is a wrapper around the :cpp:func:`lammps_extract_atom`
+    This is a wrapper around the :cpp:func:`lammps_extract_atom` function of the
-    function of the C-library interface. Its documentation includes a
+    C-library interface. Its documentation includes a list of the supported
-    list of the supported keywords and their data types.
+    keywords and their data types.  Since Python needs to know the data type to
-    Since Python needs to know the data type to be able to interpret
+    be able to interpret the result, by default, this function will try to
-    the result, by default, this function will try to auto-detect the data type
+    auto-detect the data type by asking the library. You can also force a
-    by asking the library. You can also force a specific data type by setting ``dtype``
+    specific data type by setting ``dtype`` to one of the :ref:`data type
-    to one of the :ref:`data type <py_datatype_constants>` constants defined in the
+    <py_datatype_constants>` constants defined in the :py:mod:`lammps` module.
-    :py:mod:`lammps` module.
+    This function returns ``None`` if either the keyword is not recognized, or
-    This function returns ``None`` if either the keyword is not
+    an invalid data type constant is used.
-    recognized, or an invalid data type constant is used.
 
     .. note::
 
-       While the returned arrays of per-atom data are dimensioned
+       While the returned vectors or arrays of per-atom data are dimensioned for
-       for the range [0:nmax] - as is the underlying storage -
+       the range [0:nmax] - as is the underlying storage - the data is usually
-       the data is usually only valid for the range of [0:nlocal],
+       only valid for the range of [0:nlocal], unless the property of interest
-       unless the property of interest is also updated for ghost
+       is also updated for ghost atoms.  In some cases, this depends on a LAMMPS
-       atoms.  In some cases, this depends on a LAMMPS setting, see
+       setting, see for example :doc:`comm_modify vel yes <comm_modify>`.
-       for example :doc:`comm_modify vel yes <comm_modify>`.
+       The actual size can be determined by calling
+       py:meth:`extract_atom_size() <lammps.lammps.extract_atom_size>`.
 
     :param name: name of the property
     :type name:  string
@@ -1105,6 +1135,7 @@ class lammps(object):
             ctypes.POINTER(ctypes.c_int64), ctypes.POINTER(ctypes.POINTER(ctypes.c_int64)),
             ctypes.POINTER(ctypes.c_double), ctypes.POINTER(ctypes.POINTER(ctypes.c_double)),
             or NoneType
+
     """
     if dtype == LAMMPS_AUTODETECT:
       dtype = self.extract_atom_datatype(name)
@@ -2522,3 +2553,8 @@ class lammps(object):
     newcomputeid = computeid.encode()
     idx = self.lib.lammps_find_compute_neighlist(self.lmp, newcomputeid, reqid)
     return idx
+
+# Local Variables:
+# fill-column: 80
+# End:
+  

python/lammps/numpy_wrapper.py

@@ -54,7 +54,8 @@ class numpy_wrapper:
 
   # -------------------------------------------------------------------------
 
-  def extract_atom(self, name, dtype=LAMMPS_AUTODETECT, nelem=LAMMPS_AUTODETECT, dim=LAMMPS_AUTODETECT):
+  def extract_atom(self, name, dtype=LAMMPS_AUTODETECT, nelem=LAMMPS_AUTODETECT,
+                   dim=LAMMPS_AUTODETECT):
     """Retrieve per-atom properties from LAMMPS as NumPy arrays
 
     This is a wrapper around the :py:meth:`lammps.extract_atom()` method.
@@ -63,16 +64,16 @@ class numpy_wrapper:
 
     .. note::
 
-       The returned arrays of per-atom data are by default dimensioned
+       The returned vectors or arrays of per-atom data are dimensioned
-       for the range [0:nlocal] since that data is *always* valid.  The
+       according to the return value of :py:meth:`lammps.extract_atom_size()`.
-       underlying storage for the data, however, is typically allocated
+       Except for the "mass" property, the underlying storage will always be
-       for the range of [0:nmax]. Whether there is valid data in the range
+       dimensioned for the range [0:nmax].  The actual usable data may be
-       [nlocal:nlocal+nghost] depends on whether the property of interest
+       only in the range [0:nlocal] or [0:nlocal][0:dim].  Whether there is
-       is also updated for ghost atoms.  This is not often the case.  In
+       valid data in the range [nlocal:nlocal+nghost] or [nlocal:local+nghost][0:dim]
-       some cases, it depends on a LAMMPS setting, see for example
+       depends on whether the property of interest is also updated for ghost atoms.
-       :doc:`comm_modify vel yes <comm_modify>`.  By using the optional
+       Also the value of *dim* depends on the value of *name*.  By using the optional
-       *nelem* parameter the size of the returned NumPy can be overridden.
+       *nelem* and *dim* parameters the dimensions of the returned NumPy array can
-       There is no check whether the number of elements chosen is valid.
+       be overridden.  There is no check whether the number of elements chosen is valid.
 
     :param name: name of the property
     :type name:  string
@@ -89,21 +90,10 @@ class numpy_wrapper:
       dtype = self.lmp.extract_atom_datatype(name)
 
     if nelem == LAMMPS_AUTODETECT:
-      if name == "mass":
+      nelem = self.lmp.extract_atom_size(name, LMP_SIZE_ROWS)
-        nelem = self.lmp.extract_global("ntypes") + 1
-      else:
-        nelem = self.lmp.extract_global("nlocal")
     if dim == LAMMPS_AUTODETECT:
       if dtype in (LAMMPS_INT_2D, LAMMPS_DOUBLE_2D, LAMMPS_INT64_2D):
-        # TODO add other fields
+        dim = self.lmp.extract_atom_size(name, LMP_SIZE_COLS)
-        if name in ("x", "v", "f", "x0","omega", "angmom", "torque", "csforce", "vforce", "vest"):
-          dim = 3
-        elif name == "smd_data_9":
-          dim = 9
-        elif name == "smd_stress":
-          dim = 6
-        else:
-          dim = 2
       else:
         dim = 1
 
@@ -119,37 +109,6 @@ class numpy_wrapper:
 
   # -------------------------------------------------------------------------
 
-  def extract_atom_iarray(self, name, nelem, dim=1):
-    warnings.warn("deprecated, use extract_atom instead", DeprecationWarning)
-
-    if name in ['id', 'molecule']:
-      c_int_type = self.lmp.c_tagint
-    elif name in ['image']:
-      c_int_type = self.lmp.c_imageint
-    else:
-      c_int_type = c_int
-
-    if dim == 1:
-      raw_ptr = self.lmp.extract_atom(name, LAMMPS_INT)
-    else:
-      raw_ptr = self.lmp.extract_atom(name, LAMMPS_INT_2D)
-
-    return self.iarray(c_int_type, raw_ptr, nelem, dim)
-
-  # -------------------------------------------------------------------------
-
-  def extract_atom_darray(self, name, nelem, dim=1):
-    warnings.warn("deprecated, use extract_atom instead", DeprecationWarning)
-
-    if dim == 1:
-      raw_ptr = self.lmp.extract_atom(name, LAMMPS_DOUBLE)
-    else:
-      raw_ptr = self.lmp.extract_atom(name, LAMMPS_DOUBLE_2D)
-
-    return self.darray(raw_ptr, nelem, dim)
-
-  # -------------------------------------------------------------------------
-
   def extract_compute(self, cid, cstyle, ctype):
     """Retrieve data from a LAMMPS compute
 

src/atom.cpp

@@ -2999,11 +2999,13 @@ length of the data area, and a short description.
      - N double values defined by fix property/atom array name
 
 *See also*
-   :cpp:func:`lammps_extract_atom`
+   :cpp:func:`lammps_extract_atom`,
+   :cpp:func:`lammps_extract_atom_datatype`,
+   :cpp:func:`lammps_extract_atom_size`
 
 \endverbatim
  *
- * \sa extract_datatype
+ * \sa extract_datatype, extract_size
  *
  * \param  name  string with the keyword of the desired property.
                  Typically the name of the pointer variable returned
@@ -3142,7 +3144,7 @@ void *Atom::extract(const char *name)
 
 \endverbatim
  *
- * \sa extract
+ * \sa extract extract_size
  *
  * \param  name  string with the keyword of the desired property.
  * \return       data type constant for desired property or -1 */
@@ -3177,10 +3179,14 @@ int Atom::extract_datatype(const char *name)
   if (strcmp(name,"temperature") == 0) return LAMMPS_DOUBLE;
   if (strcmp(name,"heatflow") == 0) return LAMMPS_DOUBLE;
 
+  // PERI package (and in part MACHDYN)
+
   if (strcmp(name,"vfrac") == 0) return LAMMPS_DOUBLE;
   if (strcmp(name,"s0") == 0) return LAMMPS_DOUBLE;
   if (strcmp(name,"x0") == 0) return LAMMPS_DOUBLE_2D;
 
+  // AWPMD package (and in part EFF and ELECTRODE)
+
   if (strcmp(name,"espin") == 0) return LAMMPS_INT;
   if (strcmp(name,"spin") == 0) return LAMMPS_INT;   // backwards compatibility
   if (strcmp(name,"eradius") == 0) return LAMMPS_DOUBLE;
@@ -3261,6 +3267,245 @@ int Atom::extract_datatype(const char *name)
   return -1;
 }
 
+/** Provide vector or array size info of internal data of the Atom class
+ *
+\verbatim embed:rst
+
+.. versionadded:: TBD
+
+\endverbatim
+ *
+ * \sa extract extract_datatype
+ *
+ * \param  name  string with the keyword of the desired property.
+ * \param  type  either LMP_SIZE_ROWS or LMP_SIZE_COLS for per-atom array or ignored
+ * \return       size of the vector or size of the array for the requested dimension or -1 */
+
+int Atom::extract_size(const char *name, int type)
+{
+  // --------------------------------------------------------------------
+  // 6th customization section: customize by adding new variable name
+
+  const auto datatype = extract_datatype(name);
+  const auto nall = nlocal + nghost;
+  const auto ghost_vel = comm->ghost_velocity;
+
+  if ((datatype == LAMMPS_DOUBLE_2D) || (datatype == LAMMPS_INT_2D)) {
+    if (type == LMP_SIZE_ROWS) {
+      if (strcmp(name,"x") == 0) return nall;
+      if (strcmp(name,"v") == 0) {
+        if (ghost_vel) return nall;
+        else return nlocal;
+      }
+      if (strcmp(name,"f") == 0) return nall;
+      if (strcmp(name,"mu") == 0) return nall;
+      if (strcmp(name,"omega") == 0) {
+        if (ghost_vel) return nall;
+        else return nlocal;
+      }
+      if (strcmp(name,"angmom") == 0) {
+        if (ghost_vel) return nall;
+        else return nlocal;
+      }
+      if (strcmp(name,"torque") == 0) return nlocal;
+      if (strcmp(name,"quat") == 0) {
+        if (ghost_vel) return nall;
+        else return nlocal;
+      }
+
+      // PERI package
+
+      if (strcmp(name,"x0") == 0) return nall;
+
+      // SPIN package
+
+      if (strcmp(name,"sp") == 0) return nall;
+      if (strcmp(name,"fm") == 0) return nlocal;
+      if (strcmp(name,"fm_long") == 0) return nlocal;
+
+      // AWPMD package
+
+      if (strcmp(name,"cs") == 0) {
+        if (ghost_vel) return nall;
+        else return nlocal;
+      }
+      if (strcmp(name,"csforce") == 0) return nlocal;
+      if (strcmp(name,"vforce") == 0) return nlocal;
+
+      // SPH package
+
+      if (strcmp(name,"vest") == 0) return nall;
+
+      // MACHDYN package
+
+      if (strcmp(name, "smd_data_9") == 0) return LAMMPS_DOUBLE_2D;
+      if (strcmp(name, "smd_stress") == 0) return LAMMPS_DOUBLE_2D;
+
+    } else if (type == LMP_SIZE_COLS) {
+
+      if (strcmp(name,"x") == 0) return 3;
+      if (strcmp(name,"v") == 0) return 3;
+      if (strcmp(name,"f") == 0) return 3;
+      if (strcmp(name,"mu") == 0) return 4;
+      if (strcmp(name,"omega") == 0) return 3;
+      if (strcmp(name,"angmom") == 0) return 3;
+      if (strcmp(name,"torque") == 0) return 3;
+      if (strcmp(name,"quat") == 0) return 4;
+
+      // PERI package
+
+      if (strcmp(name,"x0") == 0) return 3;
+
+      // SPIN package
+
+      if (strcmp(name,"sp") == 0) return 4;
+      if (strcmp(name,"fm") == 0) return 3;
+      if (strcmp(name,"fm_long") == 0) return 3;
+
+      // AWPMD package
+
+      if (strcmp(name,"cs") == 0) return 2;
+      if (strcmp(name,"csforce") == 0) return 2;
+      if (strcmp(name,"vforce") == 0) return 3;
+
+      // SPH package
+
+      if (strcmp(name,"vest") == 0) return 3;
+
+      // MACHDYN package
+
+      if (strcmp(name, "smd_data_9") == 0) return 9;
+      if (strcmp(name, "smd_stress") == 0) return 6;
+    }
+
+    // custom arrays
+
+    if (utils::strmatch(name,"^[id]2_")) {
+      int which = 0;
+      if (name[0] == 'd') which = 1;
+
+      int index,flag,cols,ghost;
+      index = find_custom_ghost(&name[3],flag,cols,ghost);
+
+      // consistency checks
+      if (index < 0) return -1;
+      if (which != flag) return -1;
+      if (!cols) return -1;
+
+      if (type == LMP_SIZE_ROWS) {
+        if (ghost) return nall;
+        else return nlocal;
+      } else if (type == LMP_SIZE_COLS) {
+        return cols;
+      }
+    }
+  } else {
+
+    if (strcmp(name,"mass") == 0) return ntypes + 1;
+
+    if (strcmp(name,"id") == 0) return nall;
+    if (strcmp(name,"type") == 0) return nall;
+    if (strcmp(name,"mask") == 0) return nall;
+    if (strcmp(name,"image") == 0) return nlocal;
+    if (strcmp(name,"molecule") == 0) return nall;
+    if (strcmp(name,"q") == 0) return nall;
+    if (strcmp(name,"radius") == 0) return nall;
+    if (strcmp(name,"rmass") == 0) return nall;
+
+    // ASPHERE package
+
+    if (strcmp(name,"ellipsoid") == 0) return nlocal;
+
+    // BODY package
+
+    if (strcmp(name,"line") == 0) return nlocal;
+    if (strcmp(name,"tri") == 0) return nlocal;
+    if (strcmp(name,"body") == 0) return nlocal;
+
+    // PERI package (and in part MACHDYN)
+
+    if (strcmp(name,"vfrac") == 0) return nall;
+    if (strcmp(name,"s0") == 0) return nall;
+
+    // AWPMD package (and in part EFF and ELECTRODE)
+
+    if (strcmp(name,"espin") == 0) return nall;
+    if (strcmp(name,"spin") == 0) return nall;   // backwards compatibility
+    if (strcmp(name,"eradius") == 0) return nall;
+    if (strcmp(name,"ervel") == 0) return nlocal;
+    if (strcmp(name,"erforce") == 0) return nlocal;
+    if (strcmp(name,"ervelforce") == 0) return nlocal;
+    if (strcmp(name,"etag") == 0) return nall;
+
+    // CG-DNA package
+
+    if (strcmp(name,"id5p") == 0) return nall;
+
+    // RHEO package
+
+    if (strcmp(name,"temperature") == 0) return nlocal;
+    if (strcmp(name,"heatflow") == 0) return nlocal;
+    if (strcmp(name,"rheo_status") == 0) return nall;
+    if (strcmp(name,"conductivity") == 0) return nlocal;
+    if (strcmp(name,"pressure") == 0) return nlocal;
+    if (strcmp(name,"viscosity") == 0) return nlocal;
+
+    // SPH package
+
+    if (strcmp(name,"rho") == 0) return nall;
+    if (strcmp(name,"drho") == 0) return nlocal;
+    if (strcmp(name,"esph") == 0) return nall;
+    if (strcmp(name,"desph") == 0) return nlocal;
+    if (strcmp(name,"cv") == 0) return nall;
+
+    // MACHDYN package
+
+    if (strcmp(name, "contact_radius") == 0) return nall;
+    if (strcmp(name, "eff_plastic_strain") == 0) return nlocal;
+    if (strcmp(name, "eff_plastic_strain_rate") == 0) return nlocal;
+    if (strcmp(name, "damage") == 0) return nlocal;
+
+    // DPD-REACT package
+
+    if (strcmp(name,"dpdTheta") == 0) return nall;
+
+    // DPD-MESO package
+
+    if (strcmp(name,"edpd_temp") == 0) return nall;
+
+    // DIELECTRIC package
+
+    if (strcmp(name,"area") == 0) return nall;
+    if (strcmp(name,"ed") == 0) return nall;
+    if (strcmp(name,"em") == 0) return nall;
+    if (strcmp(name,"epsilon") == 0) return nall;
+    if (strcmp(name,"curvature") == 0) return nall;
+    if (strcmp(name,"q_unscaled") == 0) return nall;
+
+    // end of customization section
+    // --------------------------------------------------------------------
+
+    // custom vectors
+
+    if (utils::strmatch(name,"^[id]_")) {
+      int which = 0;
+      if (name[0] == 'd') which = 1;
+
+      int index,flag,cols,ghost;
+      index = find_custom_ghost(&name[2],flag,cols,ghost);
+
+      // consistency checks
+      if (index < 0) return -1;
+      if (which != flag) return -1;
+      if (cols) return -1;
+
+      if (ghost) return nall;
+      else return nlocal;
+    }
+  }
+  return -1;
+}
+
 /* ----------------------------------------------------------------------
    return # of bytes of allocated memory
    call to avec tallies per-atom vectors

src/atom.h

@@ -378,6 +378,7 @@ class Atom : protected Pointers {
 
   void *extract(const char *);
   int extract_datatype(const char *);
+  int extract_size(const char *, int);
 
   inline int *get_map_array() { return map_array; };
   inline int get_map_size() { return map_tag_max + 1; };

src/library.cpp

@@ -2087,10 +2087,13 @@ int lammps_map_atom(void *handle, const void *id)
 
 .. versionadded:: 18Sep2020
 
-This function returns an integer that encodes the data type of the per-atom
+This function returns an integer that encodes the data type of the
-property with the specified name. See :cpp:enum:`_LMP_DATATYPE_CONST` for valid
+per-atom property with the specified name. See
-values. Callers of :cpp:func:`lammps_extract_atom` can use this information
+:cpp:enum:`_LMP_DATATYPE_CONST` for valid values. Callers of
-to then decide how to cast the ``void *`` pointer and access the data.
+:cpp:func:`lammps_extract_atom` can use this information to decide how
+to cast the ``void *`` pointer and access the data.  In addition,
+:cpp:func:`lammps_extract_atom_size` can be used to get information
+about the vector or array dimensions.
 
 \endverbatim
  *
@@ -2108,18 +2111,53 @@ int lammps_extract_atom_datatype(void *handle, const char *name)
 
 /* ---------------------------------------------------------------------- */
 
+/** Get dimension info of a LAMMPS per-atom property
+ *
+\verbatim embed:rst
+
+.. versionadded:: TBD
+
+This function returns an integer with the size of the per-atom
+property with the specified name.  This allows to accurately determine
+the size of the per-atom data vectors or arrays.  For per-atom arrays,
+the *type* argument is required to return either the number of rows or the
+number of columns.  It is ignored for per-atom vectors.
+
+Callers of :cpp:func:`lammps_extract_atom` can use this information in
+combination with the result from :cpp:func:`lammps_extract_atom_datatype`
+to decide how to cast the ``void *`` pointer and access the data.
+
+\endverbatim
+ *
+ * \param  handle  pointer to a previously created LAMMPS instance
+ * \param  name    string with the name of the extracted property
+ * \param  type    either LMP_SIZE_ROWS or LMP_SIZE_COLS if *name* refers
+                   to a per-atom array otherwise ignored
+ * \return         integer with the size of the vector or array dimension or -1
+ * */
+
+int lammps_extract_atom_size(void *handle, const char *name, int type)
+{
+  auto lmp = (LAMMPS *) handle;
+  return lmp->atom->extract_size(name, type);
+}
+
+/* ---------------------------------------------------------------------- */
+
 /** Get pointer to a LAMMPS per-atom property.
  *
 \verbatim embed:rst
 
-This function returns a pointer to the location of per-atom properties
+This function returns a pointer to the location of per-atom properties (and
-(and per-atom-type properties in the case of the 'mass' keyword).
+per-atom-type properties in the case of the 'mass' keyword).  Per-atom data is
-Per-atom data is distributed across sub-domains and thus MPI ranks.  The
+distributed across sub-domains and thus MPI ranks.  The returned pointer is cast
-returned pointer is cast to ``void *`` and needs to be cast to a pointer
+to ``void *`` and needs to be cast to a pointer of data type that the entity
-of data type that the entity represents.
+represents.  You can use the functions :cpp:func:`lammps_extract_atom_datatype`
+and :cpp:func:`lammps_extract_atom_size` to determine data type, dimensions and
+sizes of the storage pointed to by the returned pointer.
 
-A table with supported keywords is included in the documentation
+A table with supported keywords is included in the documentation of the
-of the :cpp:func:`Atom::extract() <LAMMPS_NS::Atom::extract>` function.
+:cpp:func:`Atom::extract() <LAMMPS_NS::Atom::extract>` function.
 
 .. warning::
 
@@ -7027,5 +7065,5 @@ int lammps_python_api_version() {
 }
 
 // Local Variables:
-// fill-column: 72
+// fill-column: 80
 // End:

src/library.h

@@ -172,6 +172,7 @@ int lammps_map_atom(void *handle, const void *id);
  * ---------------------------------------------------------------------- */
 
 int lammps_extract_atom_datatype(void *handle, const char *name);
+int lammps_extract_atom_size(void *handle, const char *name, int type);
 void *lammps_extract_atom(void *handle, const char *name);
 
 /* ----------------------------------------------------------------------

tools/swig/lammps.i

@@ -130,6 +130,7 @@ extern void  *lammps_extract_pair(void *handle, const char *name);
 extern int    lammps_map_atom(void *handle, const void *id);
 
 extern int    lammps_extract_atom_datatype(void *handle, const char *name);
+extern int    lammps_extract_atom_size(void *handle, const char *name, int type);
 extern void  *lammps_extract_atom(void *handle, const char *name);
 
 extern void  *lammps_extract_compute(void *handle, const char *id, int, int);
@@ -319,6 +320,7 @@ extern void  *lammps_extract_pair(void *handle, const char *name);
 extern int    lammps_map_atom(void *handle, const void *id);
 
 extern int    lammps_extract_atom_datatype(void *handle, const char *name);
+extern int    lammps_extract_atom_size(void *handle, const char *name, int type);
 extern void  *lammps_extract_atom(void *handle, const char *name);
 
 extern void  *lammps_extract_compute(void *handle, const char *id, int, int);

unittest/c-library/test_library_properties.cpp

@@ -695,6 +695,7 @@ TEST_F(LibraryProperties, has_error)
 class AtomProperties : public ::testing::Test {
 protected:
     void *lmp;
+    int ntypes, nlocal, nall;
 
     AtomProperties()           = default;
     ~AtomProperties() override = default;
@@ -732,6 +733,10 @@ protected:
         lammps_command(lmp, "set atom 2 i_two[2]  3");
         lammps_command(lmp, "set atom * d_four[1] -1.3");
         lammps_command(lmp, "set atom * d_four[2]  3.5");
+        ntypes = lammps_extract_setting(lmp, "ntypes");
+        nlocal = lammps_extract_setting(lmp, "nlocal");
+        nall   = lammps_extract_setting(lmp, "nall");
+
         output = ::testing::internal::GetCapturedStdout();
         if (verbose) std::cout << output;
     }
@@ -754,10 +759,12 @@ TEST_F(AtomProperties, invalid)
 TEST_F(AtomProperties, mass)
 {
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "mass"), LAMMPS_DOUBLE);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "mass", 0), ntypes + 1);
     auto *mass = (double *)lammps_extract_atom(lmp, "mass");
     ASSERT_NE(mass, nullptr);
     ASSERT_DOUBLE_EQ(mass[1], 3.0);
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "rmass"), LAMMPS_DOUBLE);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "rmass", 0), nall);
     mass = (double *)lammps_extract_atom(lmp, "rmass");
     ASSERT_NE(mass, nullptr);
     ASSERT_DOUBLE_EQ(mass[0], 2.0);
@@ -767,6 +774,7 @@ TEST_F(AtomProperties, mass)
 TEST_F(AtomProperties, charge)
 {
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "q"), LAMMPS_DOUBLE);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "rmass", 0), nall);
     auto *charge = (double *)lammps_extract_atom(lmp, "q");
     ASSERT_NE(charge, nullptr);
     ASSERT_DOUBLE_EQ(charge[0], -1.0);
@@ -776,6 +784,7 @@ TEST_F(AtomProperties, charge)
 TEST_F(AtomProperties, molecule)
 {
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "molecule"), LAMMPS_TAGINT);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "molecule", 0), nall);
     auto *molecule = (tagint *)lammps_extract_atom(lmp, "molecule");
     ASSERT_NE(molecule, nullptr);
     ASSERT_EQ(molecule[0], 2);
@@ -785,6 +794,7 @@ TEST_F(AtomProperties, molecule)
 TEST_F(AtomProperties, id)
 {
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "id"), LAMMPS_TAGINT);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "id", 0), nall);
     auto *id = (tagint *)lammps_extract_atom(lmp, "id");
     ASSERT_NE(id, nullptr);
     ASSERT_EQ(id[0], 1);
@@ -794,6 +804,7 @@ TEST_F(AtomProperties, id)
 TEST_F(AtomProperties, type)
 {
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "type"), LAMMPS_INT);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "type", 0), nall);
     int *type = (int *)lammps_extract_atom(lmp, "type");
     ASSERT_NE(type, nullptr);
     ASSERT_EQ(type[0], 1);
@@ -803,6 +814,8 @@ TEST_F(AtomProperties, type)
 TEST_F(AtomProperties, position)
 {
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "x"), LAMMPS_DOUBLE_2D);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "x", LMP_SIZE_ROWS), nall);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "x", LMP_SIZE_COLS), 3);
     auto **x = (double **)lammps_extract_atom(lmp, "x");
     ASSERT_NE(x, nullptr);
     EXPECT_DOUBLE_EQ(x[0][0], 1.0);
@@ -816,15 +829,21 @@ TEST_F(AtomProperties, position)
 TEST_F(AtomProperties, custom)
 {
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "i_one"), LAMMPS_INT);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "i_one", 0), nlocal);
     auto *one = (int *)lammps_extract_atom(lmp, "i_one");
     ASSERT_NE(one, nullptr);
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "i2_two"), LAMMPS_INT_2D);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "i2_two", LMP_SIZE_ROWS), nlocal);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "i2_two", LMP_SIZE_COLS), 2);
     auto **two = (int **)lammps_extract_atom(lmp, "i2_two");
     ASSERT_NE(two, nullptr);
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "d_three"), LAMMPS_DOUBLE);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "d_three", 0), nlocal);
     auto *three = (double *)lammps_extract_atom(lmp, "d_three");
     ASSERT_NE(three, nullptr);
     EXPECT_EQ(lammps_extract_atom_datatype(lmp, "d2_four"), LAMMPS_DOUBLE_2D);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "d2_four", LMP_SIZE_ROWS), nlocal);
+    EXPECT_EQ(lammps_extract_atom_size(lmp, "d2_four", LMP_SIZE_COLS), 2);
     auto **four = (double **)lammps_extract_atom(lmp, "d2_four");
     ASSERT_NE(four, nullptr);
 

unittest/python/python-numpy.py

@@ -155,37 +155,6 @@ class PythonNumpy(unittest.TestCase):
         self.assertEqual(values[1,0], 1.5)
         self.assertEqual(values[1,3], 1.5)
 
-    def testExtractAtomDeprecated(self):
-        self.lmp.command("units lj")
-        self.lmp.command("atom_style atomic")
-        self.lmp.command("atom_modify map array")
-        self.lmp.command("region box block 0 2 0 2 0 2")
-        self.lmp.command("create_box 1 box")
-
-        x = [
-          1.0, 1.0, 1.0,
-          1.0, 1.0, 1.5
-        ]
-
-        types = [1, 1]
-
-        self.assertEqual(self.lmp.create_atoms(2, id=None, type=types, x=x), 2)
-        nlocal = self.lmp.extract_global("nlocal", LAMMPS_INT)
-        self.assertEqual(nlocal, 2)
-
-        ident = self.lmp.numpy.extract_atom_iarray("id", nlocal, dim=1)
-        self.assertEqual(len(ident), 2)
-
-        ntypes = self.lmp.extract_global("ntypes", LAMMPS_INT)
-        self.assertEqual(ntypes, 1)
-
-        x = self.lmp.numpy.extract_atom_darray("x", nlocal, dim=3)
-        v = self.lmp.numpy.extract_atom_darray("v", nlocal, dim=3)
-        self.assertEqual(len(x), 2)
-        self.assertTrue((x[0] == (1.0, 1.0, 1.0)).all())
-        self.assertTrue((x[1] == (1.0, 1.0, 1.5)).all())
-        self.assertEqual(len(v), 2)
-
     def testExtractAtom(self):
         self.lmp.command("units lj")
         self.lmp.command("atom_style atomic")