diff --git a/doc/src/JPG/pair_body_rounded.jpg b/doc/src/JPG/pair_body_rounded.jpg
new file mode 100644
index 0000000000000000000000000000000000000000..cd745c44e27e96585be021b3dcca342f686e15e9
Binary files /dev/null and b/doc/src/JPG/pair_body_rounded.jpg differ
diff --git a/doc/src/JPG/pair_body_rounded.png b/doc/src/JPG/pair_body_rounded.png
deleted file mode 100644
index 709105ca325764adb567f123a8797ca4f76c581e..0000000000000000000000000000000000000000
Binary files a/doc/src/JPG/pair_body_rounded.png and /dev/null differ
diff --git a/doc/src/body.txt b/doc/src/body.txt
index 0e64e6ad5bad3d77105e4509d8b34468e7bbda13..be5b3d4edfd486d31c03d3fc978883e7dfacbb51 100644
--- a/doc/src/body.txt
+++ b/doc/src/body.txt
@@ -27,9 +27,9 @@ styles supported by LAMMPS are as follows.  The name in the first
 column is used as the {bstyle} argument for the "atom_style
 body"_atom_style.html command.
 
-{nparticle} | rigid body with N sub-particles |
-{rounded/polygon} | 2d polygons with N vertices :tb(c=2,s=|)
-{rounded/polyhedron} | 3d polyhedra with N vertices, E edges and F faces  :tb(c=2,s=|)
+{nparticle} : rigid body with N sub-particles
+{rounded/polygon} : 2d polygons with N vertices
+{rounded/polyhedron} : 3d polyhedra with N vertices, E edges and F faces :tb(s=:)
 
 The body style determines what attributes are stored for each body and
 thus how they can be used to compute pairwise body/body or
@@ -180,11 +180,11 @@ The {bflag2} argument is ignored.
 
 The {rounded/polygon} body style represents body particles as a 2d
 polygon with a variable number of N vertices.  This style can only be
-used for 2d models; see the "boundary"_boundary.html command.
+used for 2d models; see the "boundary"_boundary.html command.  See the
+"pair_style body/rounded/polygon" doc page for a diagram of two
+squares with rounded circles at the vertices.  Special cases for N = 1
+(circle) and N = 2 (rod with rounded ends) can also be specified.
 
-NOTE: include a diagram of a rounded polygon body particle
-
-Special cases for N = 1 (spheres) and N = 2 (rods) are also included.
 One use of this body style is for 2d discrete element models, as
 described in "Fraige"_#Fraige.
 
@@ -220,9 +220,9 @@ vertices (x1 to zN) as 3N values (with z = 0.0 for each), followed by
 2N vertex indices corresponding to the end points of the N edges,
 followed by a single diameter value = the rounded diameter of the
 circle that surrounds each vertex. The diameter value can be different
-for each body particle. These floating-point values can be
-listed on as many lines as you wish; see the
-"read_data"_read_data.html command for more details.
+for each body particle. These floating-point values can be listed on
+as many lines as you wish; see the "read_data"_read_data.html command
+for more details.
 
 The 6 moments of inertia (ixx,iyy,izz,ixy,ixz,iyz) should be the
 values consistent with the current orientation of the rigid body
@@ -235,9 +235,10 @@ position of the particle is specified by the x,y,z values in the
 {Atoms} section of the data file.
 
 For example, the following information would specify a square particle
-whose edge length is sqrt(2) and rounded diameter is 1.0 in length unit
-and initial orientation is determined from the 6 moments of inertia
-(ixx iyy izz ixy ixz iyz):
+whose edge length is sqrt(2) and rounded diameter is 1.0.  The
+orientation of the square is aligned with the xy coordinate axes which
+is consistent with the 6 moments of inertia: ixx iyy izz ixy ixz iyz =
+1 1 4 0 0 0.
 
 3 1 27
 4
@@ -265,14 +266,14 @@ body particles with a wall.
 The {rounded/polyhedron} body style represents body particles as a 3d
 polyhedron with a variable number of N vertices, E edges and F faces.
 This style can only be used for 3d models; see the
-"boundary"_boundary.html command.
+"boundary"_boundary.html command.  See the "pair_style
+body/rounded/polygon" doc page for a diagram of a two 2d squares with
+rounded circles at the vertices.  A 3d cube with rounded spheres
+at the 8 vertices and 12 rounded edges would be similar.
 
-NOTE: include a diagram of a rounded polyhedron body particle
-
-NOTE: 2d objects can also be specified as a special case, e.g.
-for a triangle, N = 3, E = 3 and F = 1.
-
-Special cases for N = 1 (spheres) and N = 2 (rods) are also valid.
+TRUNG: What are the special cases allowed for 3d, if any?  Similar to
+the N=1 (sphere) and N=2 (rod) special cases for 2d, descibed in
+previous section.
 
 This body style is for 3d discrete element models, as described in
 "Wang"_#Wang.
@@ -316,9 +317,9 @@ moments of inertia followed by the coordinates of the N vertices (x1
 to zN) as 3N values, followed by 2N vertex indices corresponding to
 the end points of the E edges, then 4*F vertex indices defining F
 faces.  The last value is the diameter value = the rounded diameter of
-the sphere that surrounds each vertex. The diameter value can be different
-for each body particle. These floating-point values
-can be listed on as many lines as you wish; see the
+the sphere that surrounds each vertex. The diameter value can be
+different for each body particle. These floating-point values can be
+listed on as many lines as you wish; see the
 "read_data"_read_data.html command for more details.  Because the
 maxmimum vertices per face is hard-coded to be 4
 (i.e. quadrilaterals), faces with more than 4 vertices need to be
@@ -340,9 +341,10 @@ position of the particle is specified by the x,y,z values in the
 {Atoms} section of the data file.
 
 For example, the following information would specify a cubic particle
-whose edge length is 2.0 and rounded diameter is 0.5 in length unit
-and initial orientation is determined from the 6 moments of inertia
-(ixx iyy izz ixy ixz iyz):
+whose edge length is 2.0 and rounded diameter is 0.5.
+The orientation of the cube is aligned with the xyz coordinate axes
+which is consistent with the 6 moments of inertia: ixx iyy izz ixy ixz
+iyz = 0.667 0.667 0.667 0 0 0.
 
 1 3 79
 8 12 6
@@ -375,6 +377,13 @@ and initial orientation is determined from the 6 moments of inertia
 3 0 4 7
 0.5 :pre
 
+The "pair_style
+body/rounded/polhedron"_pair_body_rounded_polyhedron.html command can
+be used with this body style to compute body/body interactions.  The
+"fix wall/body/polyhedron"_fix_wall_body_polygon.html command can be
+used with this body style to compute the interaction of body particles
+with a wall.
+
 :line
 
 For output purposes via the "compute
diff --git a/doc/src/pair_body_rounded_polygon.txt b/doc/src/pair_body_rounded_polygon.txt
index d611e8ec98d50a249708f27abdf0ec3f79559c9a..588a7d6ff9986609d60966e1960838ca4add4b76 100644
--- a/doc/src/pair_body_rounded_polygon.txt
+++ b/doc/src/pair_body_rounded_polygon.txt
@@ -28,9 +28,9 @@ pair_coeff 1 1 100.0 1.0 :pre
 
 Style {body/rounded/polygon} is for use with 2d models of body
 particles of style {rounded/polygon}.  It calculates pairwise
-body/body interactions as well as interactions between body and
-point-particles.  See "Section 6.14"_Section_howto.html#howto_14 of
-the manual and the "body"_body.html doc page for more details on using
+body/body interactions as well as interactions between body and point
+particles.  See "Section 6.14"_Section_howto.html#howto_14 of the
+manual and the "body"_body.html doc page for more details on using
 body particles.
 
 This pairwise interaction between rounded polygons is described in
@@ -44,36 +44,50 @@ multiple contact points.
 
 Note that when two particles interact, the effective surface of each
 polygon particle is displaced outward from each of its vertices and
-edges by half its circle diameter.  The interaction forces and
-energies between two particles are defined with respect to the
-separation of their respective rounded surfaces, not by the separation
-of the vertices and edges themselves.
+edges by half its circle diameter (as in the diagram below of a gray
+and yellow square particle).  The interaction forces and energies
+between two particles are defined with respect to the separation of
+their respective rounded surfaces, not by the separation of the
+vertices and edges themselves.
 
 This means that the specified cutoff in the pair_style command should
 be large enough to encompass the center-to-center distance between two
-particles (at any orientation) which would produce a surface-surface
-overlap.  For example, consider two square particles with edge length
-= 1.0 and circle diameter 0.2.  The maximum distance of one polygon's
-surface from its center is not sqrt(2)/2, but (sqrt(2)+0.1)/2.  Thus
-the cutoff distance should be sqrt(2) + 0.1, since the surfaces of two
-particles that far apart could be touching.
-
-The forces between vertex-vertex, vertex-edge, vertex-face and edge-edge
+particles (at any orientation) which would produce an overlap of the
+two surfaces.  For example, consider two square particles with edge
+length = 1.0 and circle diameter 0.2.  The maximum distance of one
+polygon's surface from its center is not sqrt(2)/2, but
+(sqrt(2)+0.1)/2.  Thus the cutoff distance should be sqrt(2) + 0.1,
+since the surfaces of two particles that far apart could be touching.
+
+The forces between vertex-vertex, vertex-edge, and edge-edge overlaps
 are given by:
 
 :c,image(Eqs/pair_body_rounded.jpg)
 
 :c,image(JPG/pair_body_rounded.jpg)
 
-In "Fraige"_#Fraige, the tangential friction force between two particles
-that are in contact is modeled differently prior to gross sliding
-(i.e. static friction) and during gross-sliding (kinetic friction).
-The latter takes place when the tangential deformation exceeds
-the Coulomb frictional limit.  In the current implementation, however,
-we do not take into account frictional history, i.e. we do not keep track
-of how many time steps the two particles have been in contact
-nor calculate the tangential deformation. We assume that gross sliding
-takes place right when two particles are in contact.
+In "Fraige"_#Fraige, the tangential friction force between two
+particles that are in contact is modeled differently prior to gross
+sliding (i.e. static friction) and during gross-sliding (kinetic
+friction).  The latter takes place when the tangential deformation
+exceeds the Coulomb frictional limit.  In the current implementation,
+however, we do not take into account frictional history, i.e. we do
+not keep track of how many time steps the two particles have been in
+contact nor calculate the tangential deformation.  Instead, we assume
+that gross sliding takes place as soon as two particles are in
+contact.
+
+TRUNG: The diagram label "cohesive regions" confuses me.  Are you
+saying there is some distance d for which the force is attractive,
+i.e. the particles are cohesive?  I think when d > Ri + Rj, since Ri +
+Rj is the surface/surface overlap discussed above?  If so, then the
+discussion above about the specified cutoff is wrong?  I.e. you can
+specify a large cutoff than the surface/surface overlap to get
+cohesive interactions?  If so, this should be explained up above.
+But an additional confusion is that the specied cutoff (Rc in diagram?)
+is a single number, but depedning on the orientiation of the 2
+particles they might have a suface/surface overlap at a much
+smaller value of Ri + Rj.  So what is Rc then?
 
 The following coefficients must be defined for each pair of atom types
 via the "pair_coeff"_pair_coeff.html command as in the examples above,
@@ -82,6 +96,13 @@ or in the data file read by the "read_data"_read_data.html command:
 k_n (energy/distance^2 units)
 k_na (energy/distance^2 units) :ul
 
+Effectively, k_n and k_na are the slopes of the red lines in the plot
+above for force versus distance, for r < 0 and 0 < r < rc
+respectively.  TRUNG: is this sentence correct?
+
+TRUNG: reminder to copy any change in this file
+to the pair polyhedron file as well (and vice versa)
+
 [Mixing, shift, table, tail correction, restart, rRESPA info]:
 
 This pair style does not support the "pair_modify"_pair_modify.html
diff --git a/doc/src/pair_body_rounded_polyhedron.txt b/doc/src/pair_body_rounded_polyhedron.txt
index 52f6294d8512af30fabd22e03802a447f68a896d..621254bd723413d2cd5ed14d6b3605858d7889d1 100644
--- a/doc/src/pair_body_rounded_polyhedron.txt
+++ b/doc/src/pair_body_rounded_polyhedron.txt
@@ -33,6 +33,20 @@ point-particles.  See "Section 6.14"_Section_howto.html#howto_14 of
 the manual and the "body"_body.html doc page for more details on using
 body particles.
 
+TRUNG: I think we need a paragraph here about how body/sphere
+interactions are handled.  Does this pair style only do body/body but
+allow for a body = sphere or rod or some other degenerate case?  Or
+does this pair style allow you to model a simulation of mixed body and
+point particles, where point particles are spheroids.  If so, does
+this pair style do body/body and body/point, and you use one of the
+other granular pair styles to do point/point?  I.e. a pair hybrid
+model?  Or does everything have to be defined as bodies.  Actually
+this paragraph would make more sense in the body.txt file about how to
+create a model that includes non-body particles (spheres).  And in
+this pair style file just a couple lines about which part of the
+interactions this pair style computes.  Ditto in the pair body polygon
+file.
+
 This pairwise interaction between rounded polyhedra is described in
 "Wang"_#Wang, where a polyhedron does not have sharp corners and
 edges, but is rounded at its vertices and edges by spheres centered on
@@ -59,8 +73,8 @@ surface from its center is not sqrt(3)/2, but (sqrt(3)+0.1)/2.  Thus
 the cutoff distance should be sqrt(3) + 0.1, since the surfaces of two
 particles that far apart could be touching.
 
-The forces between vertex-vertex, vertex-edge, vertex-face, edge-edge
-and edge-face are given by:
+The forces between vertex-vertex, vertex-edge, vertex-face, edge-edge,
+and edge-face overlaps are given by:
 
 :c,image(Eqs/pair_body_rounded.jpg)
 
@@ -69,12 +83,12 @@ and edge-face are given by:
 In "Wang"_#Wang, the tangential friction force between two particles
 that are in contact is modeled differently prior to gross sliding
 (i.e. static friction) and during gross-sliding (kinetic friction).
-The latter takes place when the tangential deformation exceeds
-the Coulomb frictional limit.  In the current implementation, however,
-we do not take into account frictional history, i.e. we do not keep track
-of how many time steps the two particles have been in contact
-nor calculate the tangential deformation. We assume that gross sliding
-takes place right when two particles are in contact.
+The latter takes place when the tangential deformation exceeds the
+Coulomb frictional limit.  In the current implementation, however, we
+do not take into account frictional history, i.e. we do not keep track
+of how many time steps the two particles have been in contact nor
+calculate the tangential deformation.  Instead, we assume that gross
+sliding takes place as soon as two particles are in contact.
 
 The following coefficients must be defined for each pair of atom types
 via the "pair_coeff"_pair_coeff.html command as in the examples above,
@@ -83,6 +97,10 @@ or in the data file read by the "read_data"_read_data.html command:
 k_n (energy/distance^2 units)
 k_na (energy/distance^2 units) :ul
 
+Effectively, k_n and k_na are the slopes of the red lines in the plot
+above for force versus distance, for r < 0 and 0 < r < rc
+respectively.  TRUNG: is this sentence correct?
+
 [Mixing, shift, table, tail correction, restart, rRESPA info]:
 
 This pair style does not support the "pair_modify"_pair_modify.html