diff --git a/doc/src/fix_adapt.rst b/doc/src/fix_adapt.rst
index 29eea87e64..7eaf1081f3 100644
--- a/doc/src/fix_adapt.rst
+++ b/doc/src/fix_adapt.rst
@@ -322,27 +322,31 @@ all types from 1 to :math:`N`.  A leading asterisk means all types from
 If :doc:`bond_style hybrid <bond_hybrid>` is used, *bstyle* should be a
 sub-style name. The bond styles that currently work with fix adapt are:
 
-+-----------------------------------------------------+------------+------------+
-| :doc:`class2 <bond_class2>`                         | r0         | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`fene <bond_fene>`                             | k,r0       | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`fene/nm <bond_fene>`                          | k,r0       | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`gromos <bond_gromos>`                         | k,r0       | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`harmonic <bond_harmonic>`                     | k,r0       | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`harmonic/shift <bond_harmonic_shift>`         | k,r0,r1    | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`harmonic/shift/cut <bond_harmonic_shift_cut>` | k,r0,r1    | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`harmonic/restrain <bond_harmonic_restrain>`   | k          | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`morse <bond_morse>`                           | r0         | type bonds |
-+-----------------------------------------------------+------------+------------+
-| :doc:`nonlinear <bond_nonlinear>`                   | epsilon,r0 | type bonds |
-+-----------------------------------------------------+------------+------------+
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`class2 <bond_class2>`                         | r0                        | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`fene <bond_fene>`                             | k,r0                      | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`fene/expand <bond_fene_expand>`               | k,r0,epsilon,sigma,shift  | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`fene/nm <bond_fene>`                          | k,r0                      | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`gromos <bond_gromos>`                         | k,r0                      | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`harmonic <bond_harmonic>`                     | k,r0                      | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`harmonic/restrain <bond_harmonic_restrain>`   | k                         | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`harmonic/shift <bond_harmonic_shift>`         | k,r0,r1                   | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`harmonic/shift/cut <bond_harmonic_shift_cut>` | k,r0,r1                   | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`mm3 <bond_mm3>`                               | k,r0                      | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`morse <bond_morse>`                           | r0                        | type bonds |
++-----------------------------------------------------+---------------------------+------------+
+| :doc:`nonlinear <bond_nonlinear>`                   | epsilon,r0                | type bonds |
++-----------------------------------------------------+---------------------------+------------+
 
 ----------
 
@@ -365,13 +369,31 @@ all types from 1 to :math:`N`.  A leading asterisk means all types from
 If :doc:`angle_style hybrid <angle_hybrid>` is used, *astyle* should be a
 sub-style name. The angle styles that currently work with fix adapt are:
 
-+-----------------------------------------------+----------+-------------+
-| :doc:`harmonic <angle_harmonic>`              | k,theta0 | type angles |
-+-----------------------------------------------+----------+-------------+
-| :doc:`cosine <angle_cosine>`                  | k        | type angles |
-+-----------------------------------------------+----------+-------------+
-| :doc:`cosine/squared <angle_cosine_squared>`  | k,theta0 | type angles |
-+-----------------------------------------------+----------+-------------+
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`harmonic <angle_harmonic>`                                   | k,theta0        | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`charmm <angle_charmm>`                                       | k,theta0        | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`class2 <angle_class2>`                                       | k2,k3,k4,theta0 | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`cosine <angle_cosine>`                                       | k               | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`cosine/periodic <angle_cosine_periodic>`                     | k,b,n           | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`cosine/squared/restricted <angle_cosine_squared_restricted>` | k,theta0        | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`dipole <angle_dipole>`                                       | k,gamma0        | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`fourier <angle_fourier>`                                     | k,c0,c1,c2      | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`fourier/simple <angle_fourier_simple>`                       | k,c,n           | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`mm3 <angle_mm3>`                                             | k,theta0        | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`quartic <angle_quartic>`                                     | k2,k3,k4,theta0 | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
+| :doc:`spica <angle_spica>`                                         | k,theta0        | type angles |
++--------------------------------------------------------------------+-----------------+-------------+
 
 Note that internally, theta0 is stored in radians, so the variable
 this fix uses to reset theta0 needs to generate values in radians.
diff --git a/doc/src/read_dump.rst b/doc/src/read_dump.rst
index f912390b10..a803feb65e 100644
--- a/doc/src/read_dump.rst
+++ b/doc/src/read_dump.rst
@@ -115,10 +115,11 @@ to tell LAMMPS how many parallel files exist, via its specified
 
 The format of the dump file is selected through the *format* keyword.
 If specified, it must be the last keyword used, since all remaining
-arguments are passed on to the dump reader.  The *native* format is
-for native LAMMPS dump files, written with a :doc:`dump atom <dump>`
-or :doc:`dump custom <dump>` command.  The *xyz* format is for generic XYZ
-formatted dump files.  These formats take no additional values.
+arguments are passed on to the dump reader.  The *native* format is for
+native LAMMPS dump files, written with a :doc:`dump atom <dump>` or
+:doc:`dump custom <dump>` command.  The *xyz* format is for generic XYZ
+formatted dump files (see details below).  These formats take no
+additional values.
 
 The *molfile* format supports reading data through using the `VMD <vmd_>`_
 molfile plugin interface. This dump reader format is only available,
@@ -231,17 +232,26 @@ will then have a label corresponding to the fix-ID rather than "x" or
 labels for fields *id* and *type*\ .
 
 For dump files in *xyz* format, only the *type*, *x*, *y*, and *z*
-fields are supported.  If the atoms have an integer number as label
-(like when they are created with default setting by :doc:`dump style
-<dump>` *xyz*), that number will be used at atom type. The the atom
-label is a string, it is required that a map is created using the
-:doc:`labelmap command <labelmap>` that maps those labels to the
-(numeric) atom types that LAMMPS uses.  The xyz format dump file does
-not store atom IDs, so these are assigned consecutively to the atoms as
-they appear in the dump file, starting from 1.  Thus you should ensure
-that order of atoms is consistent from snapshot to snapshot in the XYZ
-dump file.  See the :doc:`dump_modify sort <dump_modify>` command if the
-XYZ dump file was written by LAMMPS.
+fields are supported.  There are many variants of the XYZ file format.
+LAMMPS will read the number of atoms from the first line of each frame,
+ignore the second (title) line, and then read one line for each atom in the format:
+
+.. parsed-literal::
+
+   <label> <x coordinate>  <y coordinate> <z coordinate>
+
+
+If label is an integer number as label (like with XYZ files created by
+created with default settings by :doc:`dump style <dump>` *xyz*), that
+number will be used at atom type.  If the atom label is a string, then
+it is required that a map is created using the :doc:`labelmap command
+<labelmap>`. This map needs to associate each (numeric) atom type with
+one string label and that numeric type is stored.  The xyz format dump
+file does not store atom IDs, so these are assigned consecutively to the
+atoms as they appear in the dump file, starting from 1.  Thus you should
+ensure that order of atoms is consistent from snapshot to snapshot in
+the XYZ dump file.  See the :doc:`dump_modify sort <dump_modify>`
+command if the XYZ dump file was written by LAMMPS.
 
 For dump files in *molfile* format, the *x*, *y*, *z*, *vx*, *vy*, and
 *vz* fields can be specified.  However, not all molfile formats store
diff --git a/src/CG-SPICA/angle_spica.cpp b/src/CG-SPICA/angle_spica.cpp
index e315e20f13..913428cd9b 100644
--- a/src/CG-SPICA/angle_spica.cpp
+++ b/src/CG-SPICA/angle_spica.cpp
@@ -522,3 +522,15 @@ double AngleSPICA::single(int type, int i1, int i2, int i3)
   double tk = k[type] * dtheta;
   return tk*dtheta + e13;
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleSPICA::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k") == 0) return (void *) k;
+  if (strcmp(str, "theta0") == 0) return (void *) theta0;
+  return nullptr;
+}
diff --git a/src/CG-SPICA/angle_spica.h b/src/CG-SPICA/angle_spica.h
index 539512c0e9..5e590ba7a0 100644
--- a/src/CG-SPICA/angle_spica.h
+++ b/src/CG-SPICA/angle_spica.h
@@ -37,6 +37,7 @@ class AngleSPICA : public Angle {
   void read_restart(FILE *) override;
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k, *theta0;
diff --git a/src/CLASS2/angle_class2.cpp b/src/CLASS2/angle_class2.cpp
index 118179ad91..5000f9f629 100644
--- a/src/CLASS2/angle_class2.cpp
+++ b/src/CLASS2/angle_class2.cpp
@@ -467,3 +467,17 @@ double AngleClass2::single(int type, int i1, int i2, int i3)
 
   return energy;
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleClass2::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k2") == 0) return (void *) k2;
+  if (strcmp(str, "k3") == 0) return (void *) k3;
+  if (strcmp(str, "k4") == 0) return (void *) k4;
+  if (strcmp(str, "theta0") == 0) return (void *) theta0;
+  return nullptr;
+}
diff --git a/src/CLASS2/angle_class2.h b/src/CLASS2/angle_class2.h
index f5fbd62b57..4ed6f344ae 100644
--- a/src/CLASS2/angle_class2.h
+++ b/src/CLASS2/angle_class2.h
@@ -35,6 +35,7 @@ class AngleClass2 : public Angle {
   void read_restart(FILE *) override;
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *theta0, *k2, *k3, *k4;
diff --git a/src/DIPOLE/angle_dipole.cpp b/src/DIPOLE/angle_dipole.cpp
index 6ad4a0fb4c..a66f3e1042 100644
--- a/src/DIPOLE/angle_dipole.cpp
+++ b/src/DIPOLE/angle_dipole.cpp
@@ -263,3 +263,15 @@ double AngleDipole::single(int type, int iRef, int iDip, int /*iDummy*/)
 
   return kdg * deltaGamma;    // energy
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleDipole::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k") == 0) return (void *) k;
+  if (strcmp(str, "gamma0") == 0) return (void *) gamma0;
+  return nullptr;
+}
diff --git a/src/DIPOLE/angle_dipole.h b/src/DIPOLE/angle_dipole.h
index 2e55722673..de0f958f98 100644
--- a/src/DIPOLE/angle_dipole.h
+++ b/src/DIPOLE/angle_dipole.h
@@ -36,6 +36,7 @@ class AngleDipole : public Angle {
   void read_restart(FILE *) override;
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k, *gamma0;
diff --git a/src/EXTRA-MOLECULE/angle_cosine_periodic.cpp b/src/EXTRA-MOLECULE/angle_cosine_periodic.cpp
index 0b2a6d336d..6aae8d3ff4 100644
--- a/src/EXTRA-MOLECULE/angle_cosine_periodic.cpp
+++ b/src/EXTRA-MOLECULE/angle_cosine_periodic.cpp
@@ -336,3 +336,16 @@ void AngleCosinePeriodic::born_matrix(int type, int i1, int i2, int i3, double &
   du = prefactor * sin(m_angle) / s;
   du2 = prefactor * (c * sin(m_angle) - s * cos(m_angle) * multiplicity[type]) / (s * s * s);
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleCosinePeriodic::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k") == 0) return (void *) k;
+  if (strcmp(str, "b") == 0) return (void *) b;
+  if (strcmp(str, "multiplicity") == 0) return (void *) multiplicity;
+  return nullptr;
+}
diff --git a/src/EXTRA-MOLECULE/angle_cosine_periodic.h b/src/EXTRA-MOLECULE/angle_cosine_periodic.h
index f04ed04784..f63029919e 100644
--- a/src/EXTRA-MOLECULE/angle_cosine_periodic.h
+++ b/src/EXTRA-MOLECULE/angle_cosine_periodic.h
@@ -36,6 +36,7 @@ class AngleCosinePeriodic : public Angle {
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
   void born_matrix(int type, int i1, int i2, int i3, double &du, double &du2) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k;
diff --git a/src/EXTRA-MOLECULE/angle_cosine_squared_restricted.cpp b/src/EXTRA-MOLECULE/angle_cosine_squared_restricted.cpp
index 2da31ef893..c6d78ea133 100644
--- a/src/EXTRA-MOLECULE/angle_cosine_squared_restricted.cpp
+++ b/src/EXTRA-MOLECULE/angle_cosine_squared_restricted.cpp
@@ -296,3 +296,15 @@ void AngleCosineSquaredRestricted::born_matrix(int type, int i1, int i2, int i3,
 
   du2 = 2 * k[type] * numerator / denominator;
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleCosineSquaredRestricted::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k") == 0) return (void *) k;
+  if (strcmp(str, "theta0") == 0) return (void *) theta0;
+  return nullptr;
+}
diff --git a/src/EXTRA-MOLECULE/angle_cosine_squared_restricted.h b/src/EXTRA-MOLECULE/angle_cosine_squared_restricted.h
index 674252b7d0..b38b6bc4bd 100644
--- a/src/EXTRA-MOLECULE/angle_cosine_squared_restricted.h
+++ b/src/EXTRA-MOLECULE/angle_cosine_squared_restricted.h
@@ -36,6 +36,7 @@ class AngleCosineSquaredRestricted : public Angle {
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
   void born_matrix(int type, int i1, int i2, int i3, double &du, double &du2) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k, *theta0;
diff --git a/src/EXTRA-MOLECULE/angle_fourier.cpp b/src/EXTRA-MOLECULE/angle_fourier.cpp
index da1667c06f..abcda6d036 100644
--- a/src/EXTRA-MOLECULE/angle_fourier.cpp
+++ b/src/EXTRA-MOLECULE/angle_fourier.cpp
@@ -309,3 +309,16 @@ void AngleFourier::born_matrix(int type, int i1, int i2, int i3, double &du, dou
   du = k[type] * (C1[type] + 4 * C2[type] * c);
 }
 
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleFourier::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k") == 0) return (void *) k;
+  if (strcmp(str, "C0") == 0) return (void *) C0;
+  if (strcmp(str, "C1") == 0) return (void *) C1;
+  if (strcmp(str, "C2") == 0) return (void *) C2;
+  return nullptr;
+}
diff --git a/src/EXTRA-MOLECULE/angle_fourier.h b/src/EXTRA-MOLECULE/angle_fourier.h
index 8fa5d14b26..c0e30c8e1a 100644
--- a/src/EXTRA-MOLECULE/angle_fourier.h
+++ b/src/EXTRA-MOLECULE/angle_fourier.h
@@ -36,6 +36,7 @@ class AngleFourier : public Angle {
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
   void born_matrix(int type, int i1, int i2, int i3, double &du, double &du2) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k, *C0, *C1, *C2;
diff --git a/src/EXTRA-MOLECULE/angle_fourier_simple.cpp b/src/EXTRA-MOLECULE/angle_fourier_simple.cpp
index 6de7956ffa..143a008039 100644
--- a/src/EXTRA-MOLECULE/angle_fourier_simple.cpp
+++ b/src/EXTRA-MOLECULE/angle_fourier_simple.cpp
@@ -316,3 +316,16 @@ void AngleFourierSimple::born_matrix(int type, int i1, int i2, int i3, double &d
   du2 = k[type] * C[type] * N[type] * (cos(theta) * sin(N[type] * theta)
                   - N[type] * sin(theta) * cos(N[type] * theta)) / pow(sin(theta),3);
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleFourierSimple::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k") == 0) return (void *) k;
+  if (strcmp(str, "C") == 0) return (void *) C;
+  if (strcmp(str, "N") == 0) return (void *) N;
+  return nullptr;
+}
diff --git a/src/EXTRA-MOLECULE/angle_fourier_simple.h b/src/EXTRA-MOLECULE/angle_fourier_simple.h
index 3296ba6067..d37b3a83a8 100644
--- a/src/EXTRA-MOLECULE/angle_fourier_simple.h
+++ b/src/EXTRA-MOLECULE/angle_fourier_simple.h
@@ -36,6 +36,7 @@ class AngleFourierSimple : public Angle {
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
   void born_matrix(int type, int i1, int i2, int i3, double &du, double &du2) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k, *C, *N;
diff --git a/src/EXTRA-MOLECULE/angle_quartic.cpp b/src/EXTRA-MOLECULE/angle_quartic.cpp
index aade6b4534..616c81c749 100644
--- a/src/EXTRA-MOLECULE/angle_quartic.cpp
+++ b/src/EXTRA-MOLECULE/angle_quartic.cpp
@@ -325,3 +325,17 @@ void AngleQuartic::born_matrix(int type, int i1, int i2, int i3, double &du, dou
   du2 = (2.0 * k2[type] + 6.0 * k3[type] * dtheta + 12.0 * k4[type] * dtheta2) / (s*s) -
           (2.0 * k2[type] * dtheta + 3.0 * k3[type] * dtheta2 + 4.0 * k4[type] * dtheta3) * c / (s*s*s);
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleQuartic::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k2") == 0) return (void *) k2;
+  if (strcmp(str, "k3") == 0) return (void *) k3;
+  if (strcmp(str, "k4") == 0) return (void *) k4;
+  if (strcmp(str, "theta0") == 0) return (void *) theta0;
+  return nullptr;
+}
diff --git a/src/EXTRA-MOLECULE/angle_quartic.h b/src/EXTRA-MOLECULE/angle_quartic.h
index 7de51b24d1..3ff7f6f3e4 100644
--- a/src/EXTRA-MOLECULE/angle_quartic.h
+++ b/src/EXTRA-MOLECULE/angle_quartic.h
@@ -36,6 +36,7 @@ class AngleQuartic : public Angle {
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
   void born_matrix(int type, int i1, int i2, int i3, double &du, double &du2) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k2, *k3, *k4, *theta0;
diff --git a/src/MOLECULE/angle_charmm.cpp b/src/MOLECULE/angle_charmm.cpp
index 1b66260c55..11b5abd699 100644
--- a/src/MOLECULE/angle_charmm.cpp
+++ b/src/MOLECULE/angle_charmm.cpp
@@ -309,3 +309,15 @@ double AngleCharmm::single(int type, int i1, int i2, int i3)
 
   return (tk * dtheta + rk * dr);
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleCharmm::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k") == 0) return (void *) k;
+  if (strcmp(str, "theta0") == 0) return (void *) theta0;
+  return nullptr;
+}
diff --git a/src/MOLECULE/angle_charmm.h b/src/MOLECULE/angle_charmm.h
index 2a77ac7864..cc4095e90f 100644
--- a/src/MOLECULE/angle_charmm.h
+++ b/src/MOLECULE/angle_charmm.h
@@ -35,6 +35,7 @@ class AngleCharmm : public Angle {
   void read_restart(FILE *) override;
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k, *theta0, *k_ub, *r_ub;
diff --git a/src/MOLECULE/bond_fene_expand.cpp b/src/MOLECULE/bond_fene_expand.cpp
index c7821b1826..e115596eb1 100644
--- a/src/MOLECULE/bond_fene_expand.cpp
+++ b/src/MOLECULE/bond_fene_expand.cpp
@@ -273,3 +273,18 @@ double BondFENEExpand::single(int type, double rsq, int /*i*/, int /*j*/, double
 
   return eng;
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *BondFENEExpand::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k") == 0) return (void *) k;
+  if (strcmp(str, "r0") == 0) return (void *) r0;
+  if (strcmp(str, "epsilon") == 0) return (void *) epsilon;
+  if (strcmp(str, "sigma") == 0) return (void *) sigma;
+  if (strcmp(str, "shift") == 0) return (void *) shift;
+  return nullptr;
+}
diff --git a/src/MOLECULE/bond_fene_expand.h b/src/MOLECULE/bond_fene_expand.h
index cdce710ea1..13524b0972 100644
--- a/src/MOLECULE/bond_fene_expand.h
+++ b/src/MOLECULE/bond_fene_expand.h
@@ -36,6 +36,7 @@ class BondFENEExpand : public Bond {
   void read_restart(FILE *) override;
   void write_data(FILE *) override;
   double single(int, double, int, int, double &) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *k, *r0, *epsilon, *sigma, *shift;
diff --git a/src/YAFF/angle_mm3.cpp b/src/YAFF/angle_mm3.cpp
index 3ff7df1653..920041f7e9 100644
--- a/src/YAFF/angle_mm3.cpp
+++ b/src/YAFF/angle_mm3.cpp
@@ -327,3 +327,15 @@ void AngleMM3::born_matrix(int type, int i1, int i2, int i3, double &du, double
   du = -k2[type] * df / s;
   du2 = k2[type] * (d2f - df  * c / s) / (s * s) ;
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *AngleMM3::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k2") == 0) return (void *) k2;
+  if (strcmp(str, "theta0") == 0) return (void *) theta0;
+  return nullptr;
+}
diff --git a/src/YAFF/angle_mm3.h b/src/YAFF/angle_mm3.h
index 22f5bd746c..126b275e72 100644
--- a/src/YAFF/angle_mm3.h
+++ b/src/YAFF/angle_mm3.h
@@ -36,6 +36,7 @@ class AngleMM3 : public Angle {
   void write_data(FILE *) override;
   double single(int, int, int, int) override;
   void born_matrix(int type, int i1, int i2, int i3, double &du, double &du2) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *theta0, *k2;
diff --git a/src/YAFF/bond_mm3.cpp b/src/YAFF/bond_mm3.cpp
index 31ce2dad3e..b3e69881e1 100644
--- a/src/YAFF/bond_mm3.cpp
+++ b/src/YAFF/bond_mm3.cpp
@@ -238,3 +238,15 @@ void BondMM3::born_matrix(int type, double rsq, int /*i*/, int /*j*/, double &du
   du = 2.0 * k2[type] * dr + 3.0 * K3 * dr2 + 4.0 * K4 * dr3;
   du2 = 2.0 * k2[type] + 6.0 * K3 * dr + 12.0 * K4 * dr2;
 }
+
+/* ----------------------------------------------------------------------
+   return ptr to internal members upon request
+------------------------------------------------------------------------ */
+
+void *BondMM3::extract(const char *str, int &dim)
+{
+  dim = 1;
+  if (strcmp(str, "k2") == 0) return (void *) k2;
+  if (strcmp(str, "r0") == 0) return (void *) r0;
+  return nullptr;
+}
diff --git a/src/YAFF/bond_mm3.h b/src/YAFF/bond_mm3.h
index ea89ac826d..b9ebf464bb 100644
--- a/src/YAFF/bond_mm3.h
+++ b/src/YAFF/bond_mm3.h
@@ -36,6 +36,7 @@ class BondMM3 : public Bond {
   void write_data(FILE *) override;
   double single(int, double, int, int, double &) override;
   void born_matrix(int, double, int, int, double &, double &) override;
+  void *extract(const char *, int &) override;
 
  protected:
   double *r0, *k2;
diff --git a/src/library.cpp b/src/library.cpp
index 75d6fecc18..0065f10454 100644
--- a/src/library.cpp
+++ b/src/library.cpp
@@ -509,10 +509,9 @@ void lammps_error(void *handle, int error_type, const char *error_text)
 
 This function tells LAMMPS to expand the string in *cmd* like it would process
 an input line fed to :cpp:func:`lammps_command` **without** executing it.
-The entire string is considered as command and need not have a
-(final) newline character.  Newline characters in the body of the
-string, however, will be treated as part of the command and will **not**
-start a second command.
+The *entire* string is considered as input and need not have a (final) newline
+character.  Newline characters in the body of the string, however, will be
+treated as part of the command and will **not** start a second command.
 
 The function returns the expanded string in a new string buffer that
 must be freed with :cpp:func:`lammps_free` after use to avoid a memory leak.
diff --git a/unittest/force-styles/tests/angle-charmm.yaml b/unittest/force-styles/tests/angle-charmm.yaml
index 52e78abaf6..8ff5118390 100644
--- a/unittest/force-styles/tests/angle-charmm.yaml
+++ b/unittest/force-styles/tests/angle-charmm.yaml
@@ -15,7 +15,9 @@ angle_coeff: ! |
   3  40.0 120.0  35.0 2.410
   4  33.0 108.5  30.0 2.163
 equilibrium: 4 1.9216075064457567 1.9425514574696887 2.0943951023931953 1.8936822384138476
-extract: ! ""
+extract: ! |
+  k 1
+  theta0 1
 natoms: 29
 init_energy: 85.42486388459771
 init_stress: ! |2-
diff --git a/unittest/force-styles/tests/angle-class2.yaml b/unittest/force-styles/tests/angle-class2.yaml
index ae2e3ff5ee..8901157d17 100644
--- a/unittest/force-styles/tests/angle-class2.yaml
+++ b/unittest/force-styles/tests/angle-class2.yaml
@@ -23,7 +23,11 @@ angle_coeff: ! |
   3  ba  10.0 10.0 1.5 1.5
   4  ba   0.0 20.0 1.5 1.5
 equilibrium: 4 1.9216075064457565 1.9373154697137058 2.0943951023931953 1.8936822384138474
-extract: ! ""
+extract: ! |
+  k2 1
+  k3 1
+  k4 1
+  theta0 1
 natoms: 29
 init_energy: 46.44089683774903
 init_stress: ! |2-
diff --git a/unittest/force-styles/tests/angle-cosine_periodic.yaml b/unittest/force-styles/tests/angle-cosine_periodic.yaml
index 5c8227fcbd..ee3e5c1469 100644
--- a/unittest/force-styles/tests/angle-cosine_periodic.yaml
+++ b/unittest/force-styles/tests/angle-cosine_periodic.yaml
@@ -15,7 +15,10 @@ angle_coeff: ! |
   3  50.0 -1 3
   4 100.0 -1 4
 equilibrium: 4 3.141592653589793 3.141592653589793 2.0943951023931957 2.356194490192345
-extract: ! ""
+extract: ! |
+  k 1
+  b 1
+  multiplicity 1
 natoms: 29
 init_energy: 1178.5476942873006
 init_stress: ! |2-
diff --git a/unittest/force-styles/tests/angle-cosine_squared_restricted.yaml b/unittest/force-styles/tests/angle-cosine_squared_restricted.yaml
index 400babb3c0..341ccb3919 100644
--- a/unittest/force-styles/tests/angle-cosine_squared_restricted.yaml
+++ b/unittest/force-styles/tests/angle-cosine_squared_restricted.yaml
@@ -17,7 +17,9 @@ angle_coeff: ! |
   3  50.0 120.0
   4 100.0 108.5
 equilibrium: 4 1.9216075064457567 1.9373154697137058 2.0943951023931953 1.8936822384138476
-extract: ! ""
+extract: ! |
+  k 1
+  theta0 1
 natoms: 29
 init_energy: 43.16721849625078
 init_stress: ! |2-
diff --git a/unittest/force-styles/tests/angle-dipole.yaml b/unittest/force-styles/tests/angle-dipole.yaml
index 877ffa19c7..3914973f0d 100644
--- a/unittest/force-styles/tests/angle-dipole.yaml
+++ b/unittest/force-styles/tests/angle-dipole.yaml
@@ -20,7 +20,9 @@ angle_coeff: ! |
   3  50.0 120.0
   4 100.0 108.5
 equilibrium: 4 1.9216075064457565 1.9373154697137058 2.0943951023931953 1.8936822384138474
-extract: ! ""
+extract: ! |
+  k 1
+  gamma0 1
 natoms: 29
 init_energy: 1003.6681304854917
 init_stress: ! |2-
diff --git a/unittest/force-styles/tests/angle-fourier.yaml b/unittest/force-styles/tests/angle-fourier.yaml
index 61165c5a92..275d62beca 100644
--- a/unittest/force-styles/tests/angle-fourier.yaml
+++ b/unittest/force-styles/tests/angle-fourier.yaml
@@ -16,7 +16,11 @@ angle_coeff: ! |
   3  50.0 0.0 0.0 1.0
   4 100.0 0.3 0.3 0.3
 equilibrium: 4 3.141592653589793 1.5707963267948966 1.5707963267948966 1.8234765819369754
-extract: ! ""
+extract: ! |
+  k 1
+  C0 1
+  C1 1
+  C2 1
 natoms: 29
 init_energy: 400.84036632010225
 init_stress: ! |-
diff --git a/unittest/force-styles/tests/angle-fourier_simple.yaml b/unittest/force-styles/tests/angle-fourier_simple.yaml
index e1a394ee3a..bd72e67912 100644
--- a/unittest/force-styles/tests/angle-fourier_simple.yaml
+++ b/unittest/force-styles/tests/angle-fourier_simple.yaml
@@ -16,7 +16,10 @@ angle_coeff: ! |
   3  50.0  1.0 3.0
   4 100.0 -0.5 1.5
 equilibrium: 4 3.141592653589793 1.5707963267948966 1.0471975511965976 2.0943951023931953
-extract: ! ""
+extract: ! |
+  k 1
+  C 1
+  N 1
 natoms: 29
 init_energy: 2474.0748013590646
 init_stress: ! |-
diff --git a/unittest/force-styles/tests/angle-mm3.yaml b/unittest/force-styles/tests/angle-mm3.yaml
index 9fb9460183..381f43187e 100644
--- a/unittest/force-styles/tests/angle-mm3.yaml
+++ b/unittest/force-styles/tests/angle-mm3.yaml
@@ -16,7 +16,9 @@ angle_coeff: ! |
   3  50.0 120.0
   4 100.0 108.5
 equilibrium: 4 1.9216075064457565 1.9373154697137058 2.0943951023931953 1.8936822384138474
-extract: ! ""
+extract: ! |
+  k2 1
+  theta0 1
 natoms: 29
 init_energy: 44.72461548562619
 init_stress: ! |2-
diff --git a/unittest/force-styles/tests/angle-quartic.yaml b/unittest/force-styles/tests/angle-quartic.yaml
index 6ded709e84..15ec06d82a 100644
--- a/unittest/force-styles/tests/angle-quartic.yaml
+++ b/unittest/force-styles/tests/angle-quartic.yaml
@@ -16,7 +16,11 @@ angle_coeff: ! |
   3 120.0  50.0  -9.5 -1.5
   4 108.5 100.0   5.0 -2.0
 equilibrium: 4 1.9216075064457565 1.9373154697137058 2.0943951023931953 1.8936822384138474
-extract: ! ""
+extract: ! |
+  k2 1
+  k3 1
+  k4 1
+  theta0 1
 natoms: 29
 init_energy: 41.0458477552901
 init_stress: ! |2-
diff --git a/unittest/force-styles/tests/angle-spica.yaml b/unittest/force-styles/tests/angle-spica.yaml
index 7f88553c70..46e8238349 100644
--- a/unittest/force-styles/tests/angle-spica.yaml
+++ b/unittest/force-styles/tests/angle-spica.yaml
@@ -20,7 +20,9 @@ angle_coeff: ! |
   3  40.0 120.0
   4  33.0 108.5
 equilibrium: 4 1.9216075064457565 1.9425514574696887 2.0943951023931953 1.8936822384138474
-extract: ! ""
+extract: ! |
+  k 1
+  theta0 1
 natoms: 29
 init_energy: 38.36438529349082
 init_stress: ! |2-
diff --git a/unittest/force-styles/tests/bond-fene_expand.yaml b/unittest/force-styles/tests/bond-fene_expand.yaml
index fc859d477c..250f89af15 100644
--- a/unittest/force-styles/tests/bond-fene_expand.yaml
+++ b/unittest/force-styles/tests/bond-fene_expand.yaml
@@ -17,7 +17,12 @@ bond_coeff: ! |
   4 650 2.4 0.015 1.2 0.15
   5 450 2   0.018 1   0.09
 equilibrium: 5 1.5550000000000002 1.117 1.321 1.3139999999999998 1.06
-extract: ! ""
+extract: ! |
+  k 1
+  r0 1
+  epsilon 1
+  sigma 1
+  shift 1
 natoms: 29
 init_energy: 5926.020859124294
 init_stress: ! |-
diff --git a/unittest/force-styles/tests/bond-mm3.yaml b/unittest/force-styles/tests/bond-mm3.yaml
index eb7443f1c2..f5ba5c237c 100644
--- a/unittest/force-styles/tests/bond-mm3.yaml
+++ b/unittest/force-styles/tests/bond-mm3.yaml
@@ -17,7 +17,9 @@ bond_coeff: ! |
   4 650.0 1.2
   5 450.0 1.0
 equilibrium: 5 1.5 1.1 1.3 1.2 1
-extract: ! ""
+extract: ! |
+  k2 1
+  r0 1
 natoms: 29
 init_energy: 4.247265008273143
 init_stress: ! |-