diff --git a/doc/src/compute_pressure_cylinder.txt b/doc/src/compute_pressure_cylinder.txt
new file mode 100644
index 0000000000000000000000000000000000000000..4865fe8a795f961375d5486b1ed1b2b5358d42a7
--- /dev/null
+++ b/doc/src/compute_pressure_cylinder.txt
@@ -0,0 +1,81 @@
+"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
+
+:link(lws,http://lammps.sandia.gov)
+:link(ld,Manual.html)
+:link(lc,Section_commands.html#comm)
+
+:line
+
+compute pressure/cylinder command :h3
+
+[Syntax:]
+
+compute ID group-ID pressure/cylinder zlo zhi Rmax bin_width :pre
+
+ID, group-ID are documented in "compute"_compute.html command
+pressure/cylinder = style name of this compute command
+zlo = minimum z-boundary for cylinder
+zhi = maximum z-boundary for cylinder
+Rmax = maximum radius to perform calculation to
+bin_width = width of radial bins to use for calculation :ul
+
+[Examples:]
+
+compute 1 all pressure/cylinder -10.0 10.0 15.0 0.25 :pre
+
+[Description:]
+
+Define a computation that calculates the pressure tensor of a system in
+cylindrical coordinates, as discussed in "(Addington)"_#Addington1.
+This is useful for systems with a single axis of rotational symmetry,
+such as cylindrical micelles or carbon nanotubes. The compute splits the
+system into radial, cylindrical-shell-type bins of width bin_width,
+centered at x=0,y=0, and calculates the radial (P_rhorho), azimuthal
+(P_phiphi), and axial (P_zz) components of the configurational pressure
+tensor. The local density is also calculated for each bin, so that the
+true pressure can be recovered as P_kin+P_conf=density*k*T+P_conf. The
+output is a global array with 5 columns; one each for bin radius, local
+number density, P_rhorho, P_phiphi, and P_zz. The number of rows is
+governed by the values of Rmax and bin_width. Pressure tensor values are
+output in pressure units.
+
+[Output info:]
+
+This compute calculates a global array with 5 columns and Rmax/bin_width
+rows. The output columns are: R (distance units), number density (inverse
+volume units), configurational radial pressure (pressure units),
+configurational azimuthal pressure (pressure units), and configurational
+axial pressure (pressure units).
+
+The values calculated by this compute are
+"intensive". The pressure values will be in pressure
+"units"_units.html. The number density values will be in
+inverse volume "units"_units.html.
+
+[Restrictions:]
+
+This compute currently calculates the pressure tensor contributions
+for pair styles only (i.e. no bond, angle, dihedral, etc. contributions
+and in the presence of bonded interactions, the result will be incorrect
+due to exclusions for special bonds) and requires pair-wise force
+calculations not available for most manybody pair styles. K-space
+calculations are also excluded. Note that this pressure compute outputs
+the configurational terms only; the kinetic contribution is not included
+and may be calculated from the number density output by P_kin=density*k*T.
+
+This compute is part of the USER-MISC package. It is only enabled
+if LAMMPS was built with that package. See the "Build
+package"_Build_package.html doc page for more info.
+
+[Related commands:]
+
+"compute temp"_compute_temp.html, "compute
+stress/atom"_compute_stress_atom.html,
+"thermo_style"_thermo_style.html,
+
+[Default:] none
+
+:line
+
+:link(Addington1)
+[(Addington)] Addington, Long, Gubbins, J Chem Phys, 149, 084109 (2018).
diff --git a/doc/src/computes.txt b/doc/src/computes.txt
index 2e2fde7fb42b95f40606ccc7b26e8ce58de38779..c836fb6cb1c980682373dd45654ce9fc535698b1 100644
--- a/doc/src/computes.txt
+++ b/doc/src/computes.txt
@@ -67,6 +67,7 @@ Computes :h1
compute_pe_atom
compute_plasticity_atom
compute_pressure
+ compute_pressure_cylinder
compute_pressure_uef
compute_property_atom
compute_property_chunk
diff --git a/doc/src/lammps.book b/doc/src/lammps.book
index 3ce6484bd9d0a1b285cfb3f58ce63b178a3dfa43..c006e77d65a0963f96376d0f5a901fb6830bb14b 100644
--- a/doc/src/lammps.book
+++ b/doc/src/lammps.book
@@ -465,6 +465,7 @@ compute_pe.html
compute_pe_atom.html
compute_plasticity_atom.html
compute_pressure.html
+compute_pressure_cylinder.html
compute_pressure_uef.html
compute_property_atom.html
compute_property_chunk.html
diff --git a/src/.gitignore b/src/.gitignore
index e3c273ebefcbebc6d23495c3fdff0616008b7e00..5699a87089cea119b080ce2b67a223fcf67ff45f 100644
--- a/src/.gitignore
+++ b/src/.gitignore
@@ -301,6 +301,8 @@
/compute_plasticity_atom.h
/compute_pressure_bocs.cpp
/compute_pressure_bocs.h
+/compute_pressure_cylinder.cpp
+/compute_pressure_cylinder.h
/compute_pressure_grem.cpp
/compute_pressure_grem.h
/compute_rigid_local.cpp
diff --git a/src/USER-MISC/README b/src/USER-MISC/README
index b3c3e859480b0ba3e6bfbb1a501d684f02d284b3..cf948df233cb53819d220e5286454d9e21319888 100644
--- a/src/USER-MISC/README
+++ b/src/USER-MISC/README
@@ -30,6 +30,7 @@ compute ackland/atom, Gerolf Ziegenhain, gerolf at ziegenhain.com, 4 Oct 2007
compute basal/atom, Christopher Barrett, cdb333 at cavs.msstate.edu, 3 Mar 2013
compute cnp/atom, Paulo Branicio (USC), branicio at usc.edu, 31 May 2017
compute entropy/atom, Pablo Piaggi (EPFL), pablo.piaggi at epfl.ch, 15 June 2018
+compute pressure/cylinder, Cody K. Addington (NCSU), , 2 Oct 2018
compute stress/mop, Romain Vermorel (U Pau) & Laurent Joly (U Lyon), romain.vermorel at univ-pau.fr & ljoly.ulyon at gmail.com, 5 Sep 18
compute stress/mop/profile, Romain Vermorel (U Pau) & Laurent Joly (U Lyon), romain.vermorel at univ-pau.fr & ljoly.ulyon at gmail.com, 5 Sep 18
compute temp/rotate, Laurent Joly (U Lyon), ljoly.ulyon at gmail.com, 8 Aug 11
diff --git a/src/USER-MISC/compute_pressure_cylinder.cpp b/src/USER-MISC/compute_pressure_cylinder.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..6fcd343e38d643163a2cb7b4182eddf44f73f169
--- /dev/null
+++ b/src/USER-MISC/compute_pressure_cylinder.cpp
@@ -0,0 +1,484 @@
+/* ----------------------------------------------------------------------
+ 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 <cmath>
+#include <cstring>
+#include <cstdlib>
+#include "compute_pressure_cylinder.h"
+#include "atom.h"
+#include "update.h"
+#include "force.h"
+#include "pair.h"
+#include "neighbor.h"
+#include "neigh_request.h"
+#include "neigh_list.h"
+#include "group.h"
+#include "memory.h"
+#include "error.h"
+#include "citeme.h"
+#include "domain.h"
+#include "math_const.h"
+
+using namespace LAMMPS_NS;
+using namespace MathConst;
+
+static const char cite_compute_pressure_cylinder[] =
+ "compute pressure/cylinder:\n\n"
+ "@Article{Addington,\n"
+ " author = {C. K. Addington, Y. Long, K. E. Gubbins},\n"
+ " title = {The pressure in interfaces having cylindrical geometry},\n"
+ " journal = {J.~Chem.~Phys.},\n"
+ " year = 2018,\n"
+ " volume = 149,\n"
+ " pages = {084109}\n"
+ "}\n\n";
+
+/*+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+ Calculate the configurational components of the pressure tensor in
+ cylindrical geometry, according to the formulation of Addington et al. (2018)
+ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
+
+ComputePressureCyl::ComputePressureCyl(LAMMPS *lmp, int narg, char **arg) :
+ Compute(lmp, narg, arg),
+ Pr_temp(NULL), Pr_all(NULL), Pz_temp(NULL), Pz_all(NULL), Pphi_temp(NULL),
+ Pphi_all(NULL), R(NULL), Rinv(NULL), R2(NULL), PrAinv(NULL), PzAinv(NULL),
+ R2kin(NULL), density_temp(NULL), invVbin(NULL), density_all(NULL),
+ tangent(NULL), ephi_x(NULL), ephi_y(NULL), binz(NULL)
+{
+ if (lmp->citeme) lmp->citeme->add(cite_compute_pressure_cylinder);
+ if (narg != 7) error->all(FLERR,"Illegal compute pressure/cylinder command");
+
+ zlo=force->numeric(FLERR,arg[3]);
+ zhi=force->numeric(FLERR,arg[4]);
+ Rmax=force->numeric(FLERR,arg[5]);
+ bin_width=force->numeric(FLERR,arg[6]);
+
+ if ((bin_width <= 0.0) || (bin_width < Rmax))
+ error->all(FLERR,"Illegal compute pressure/cylinder command");
+ if ((zhi < zlo) || ((zhi-zlo) < bin_width))
+ error->all(FLERR,"Illegal compute pressure/cylinder command");
+ if ((zhi > domain->boxhi[2]) || (zlo < domain->boxlo[2]))
+ error->all(FLERR,"Illegal compute pressure/cylinder command");
+
+ nbins=(int)(Rmax/bin_width);
+ nzbins=(int)((zhi-zlo)/bin_width);
+
+ // NOTE: at 2^22 = 4.2M bins, we will be close to exhausting allocatable
+ // memory on a 32-bit environment. so we use this as an upper limit.
+
+ if ((nbins < 1) || (nzbins < 1) || (nbins > 2>>22) || (nbins > 2>>22))
+ error->all(FLERR,"Illegal compute pressure/cylinder command");
+
+ array_flag=1;
+ vector_flag=0;
+ extarray=0;
+ size_array_cols = 5; // r, number density, Pr, Pphi, Pz
+ size_array_rows = nbins;
+
+ Pr_temp = new double[nbins];
+ Pr_all = new double[nbins];
+ Pz_temp = new double[nbins];
+ Pz_all = new double[nbins];
+ Pphi_temp = new double[nbins];
+ Pphi_all = new double[nbins];
+ R = new double[nbins];
+ R2 = new double[nbins];
+ PrAinv = new double[nbins];
+ PzAinv = new double[nbins];
+ Rinv = new double[nbins];
+ binz = new double[nzbins];
+
+ R2kin = new double[nbins];
+ density_temp = new double[nbins];
+ invVbin = new double[nbins];
+ density_all = new double[nbins];
+
+ memory->create(array,nbins,5,"PN:array");
+
+ nphi=360;
+ tangent = new double[nphi];
+ ephi_x = new double[nphi];
+ ephi_y = new double[nphi];
+
+ nktv2p = force->nktv2p;
+
+}
+
+/* ---------------------------------------------------------------------- */
+
+ComputePressureCyl::~ComputePressureCyl()
+{
+ // count all of these for memory usage
+ memory->destroy(array);
+ delete [] R;
+ delete [] Rinv;
+ delete [] R2;
+ delete [] R2kin;
+ delete [] invVbin;
+ delete [] density_temp;
+ delete [] density_all;
+ delete [] tangent;
+ delete [] ephi_x;
+ delete [] ephi_y;
+ delete [] Pr_temp;
+ delete [] Pr_all;
+ delete [] Pz_temp;
+ delete [] Pz_all;
+ delete [] Pphi_temp;
+ delete [] Pphi_all;
+ delete [] PrAinv;
+ delete [] PzAinv;
+ delete [] binz;
+}
+
+/* ---------------------------------------------------------------------- */
+
+void ComputePressureCyl::init()
+{
+ if (force->pair == NULL)
+ error->all(FLERR,"No pair style is defined for compute pressure/cylinder");
+ if (force->pair->single_enable == 0)
+ error->all(FLERR,"Pair style does not support compute pressure/cylinder");
+
+ double phi;
+
+ for (int iphi = 0; iphi < nphi; iphi++) {
+ phi=((double)iphi)*MY_PI/180.0;
+ tangent[iphi]=tan(phi);
+ ephi_x[iphi]=-sin(phi);
+ ephi_y[iphi]=cos(phi);
+ }
+
+ for (int iq = 0; iq < nbins; iq++) {
+ R[iq]=((double)iq+0.5)*bin_width;
+ Rinv[iq]=1.0/R[iq];
+ R2[iq]=R[iq]*R[iq];
+ R2kin[iq]=(((double)iq)+1.0)*bin_width;
+ R2kin[iq]*=R2kin[iq];
+ PrAinv[iq]=1.0/(2.0*MY_PI*(zhi-zlo)*R[iq]);
+ }
+ PphiAinv=1.0/((zhi-zlo)*bin_width*2.0*(double)nphi);
+
+ invVbin[0]=1.0/((zhi-zlo)*MY_PI*R2kin[0]);
+ PzAinv[0]=1.0/(MY_PI*R2kin[0]*((double)nzbins));
+
+ for (int jq = 1; jq < nbins; jq++) {
+ invVbin[jq]=1.0/((zhi-zlo)*MY_PI*(R2kin[jq]-R2kin[jq-1]));
+ PzAinv[jq]=1.0/(MY_PI*(R2kin[jq]-R2kin[jq-1])*((double)nzbins));
+ }
+
+ // need an occasional half neighbor list
+ int irequest = neighbor->request(this,instance_me);
+ neighbor->requests[irequest]->pair = 0;
+ neighbor->requests[irequest]->compute = 1;
+ neighbor->requests[irequest]->occasional = 1;
+
+ for (int zzz = 0; zzz < nzbins; zzz++) binz[zzz]=(((double)zzz)+0.5)*bin_width+zlo;
+
+}
+
+/* ---------------------------------------------------------------------- */
+
+void ComputePressureCyl::init_list(int /* id */, NeighList *ptr)
+{
+ list = ptr;
+}
+
+/* ---------------------------------------------------------------------- */
+
+
+/* ----------------------------------------------------------------------
+ count pairs and compute pair info on this proc
+ only count pair once if newton_pair is off
+ both atom I,J must be in group
+ if flag is set, compute requested info about pair
+------------------------------------------------------------------------- */
+
+void ComputePressureCyl::compute_array()
+{
+ invoked_array = update->ntimestep;
+
+ int ibin;
+
+ // clear pressures
+ for (ibin = 0; ibin < nbins; ibin++) {
+ density_temp[ibin]=0.0;
+ density_all[ibin]=0.0;
+ Pr_temp[ibin]=0.0;
+ Pr_all[ibin]=0.0;
+ Pphi_temp[ibin]=0.0;
+ Pphi_all[ibin]=0.0;
+ Pz_temp[ibin]=0.0;
+ Pz_all[ibin]=0.0;
+ }
+
+ // what processor am I?
+ int me;
+ MPI_Comm_rank(world,&me);
+
+ int i,j,ii,jj,inum,jnum,itype,jtype;
+ tagint itag,jtag;
+ double xtmp,ytmp,ztmp,delx,dely,delz;
+ double rsq,fpair,factor_coul,factor_lj;
+ int *ilist,*jlist,*numneigh,**firstneigh;
+
+ double **x = atom->x;
+ tagint *tag = atom->tag;
+ int *type = atom->type;
+ int *mask = atom->mask;
+ int nlocal = atom->nlocal;
+ double *special_coul = force->special_coul;
+ double *special_lj = force->special_lj;
+ int newton_pair = force->newton_pair;
+
+ // invoke half neighbor list (will copy or build if necessary)
+ neighbor->build_one(list);
+
+ inum = list->inum;
+ ilist = list->ilist;
+ numneigh = list->numneigh;
+ firstneigh = list->firstneigh;
+
+ // calculate number density (by radius)
+ double temp_R2;
+ for (i = 0; i < nlocal; i++) if ((x[i][2] < zhi) && (x[i][2] > zlo)) {
+ temp_R2=x[i][0]*x[i][0]+x[i][1]*x[i][1];
+ if (temp_R2 > R2kin[nbins-1]) continue; // outside of Rmax
+
+ for (j = 0; j < nbins; j++) if (temp_R2 < R2kin[j]) break;
+
+ density_temp[j]+=invVbin[j];
+ }
+ MPI_Allreduce(density_temp,density_all,nbins,MPI_DOUBLE,MPI_SUM,world);
+ for (i = 0; i < nbins; i++) array[i][1]=density_all[i]; // NEW
+
+ // loop over neighbors of my atoms
+ // skip if I or J are not in group
+ // for newton = 0 and J = ghost atom,
+ // need to insure I,J pair is only output by one proc
+ // use same itag,jtag logic as in Neighbor::neigh_half_nsq()
+ // for flag = 0, just count pair interactions within force cutoff
+ // for flag = 1, calculate requested output fields
+
+ Pair *pair = force->pair;
+ double **cutsq = force->pair->cutsq;
+
+ double r1=0.0;
+ double r2=0.0;
+ double risq,rjsq;
+ double A,B,C,D;
+ double alpha1,alpha2;
+ double xi,yi,zi,dx,dy,dz;
+ double xR,yR,zR,fn;
+ double alpha,xL,yL,zL,L2,ftphi,ftz;
+ double sqrtD;
+
+ for (ii = 0; ii < inum; ii++) {
+ i = ilist[ii];
+ if (!(mask[i] & groupbit)) continue;
+
+ xtmp = x[i][0];
+ ytmp = x[i][1];
+ ztmp = x[i][2];
+ itag = tag[i];
+ itype = type[i];
+ jlist = firstneigh[i];
+ jnum = numneigh[i];
+
+ r1=x[i][0]*x[i][0]+x[i][1]*x[i][1];
+
+ for (jj = 0; jj < jnum; jj++) {
+ j = jlist[jj];
+ factor_lj = special_lj[sbmask(j)];
+ factor_coul = special_coul[sbmask(j)];
+ j &= NEIGHMASK;
+
+ if (!(mask[j] & groupbit)) continue;
+
+ // itag = jtag is possible for long cutoffs that include images of self
+ // do calculation only on appropriate processor
+ if (newton_pair == 0 && j >= nlocal) {
+ jtag = tag[j];
+ if (itag > jtag) {
+ if ((itag+jtag) % 2 == 0) continue;
+ } else if (itag < jtag) {
+ if ((itag+jtag) % 2 == 1) continue;
+ } else {
+ if (x[j][2] < ztmp) continue;
+ if (x[j][2] == ztmp) {
+ if (x[j][1] < ytmp) continue;
+ if (x[j][1] == ytmp && x[j][0] < xtmp) continue;
+ }
+ }
+ }
+
+ delx = xtmp - x[j][0];
+ dely = ytmp - x[j][1];
+ delz = ztmp - x[j][2];
+
+ r2=x[j][0]*x[j][0]+x[j][1]*x[j][1];
+
+ // ri is smaller of r1 and r2
+ if (r2 < r1) {
+ risq=r2;
+ rjsq=r1;
+ xi=x[j][0];
+ yi=x[j][1];
+ zi=x[j][2];
+ dx=x[i][0]-x[j][0];
+ dy=x[i][1]-x[j][1];
+ dz=x[i][2]-x[j][2];
+ } else {
+ risq=r1;
+ rjsq=r2;
+ xi=x[i][0];
+ yi=x[i][1];
+ zi=x[i][2];
+ dx=x[j][0]-x[i][0];
+ dy=x[j][1]-x[i][1];
+ dz=x[j][2]-x[i][2];
+ }
+
+ rsq = delx*delx + dely*dely + delz*delz;
+ jtype = type[j];
+ if (rsq >= cutsq[itype][jtype]) continue;
+
+ pair->single(i,j,itype,jtype,rsq,factor_coul,factor_lj,fpair);
+
+ A=dx*dx+dy*dy;
+ B=2.0*(xi*dx+yi*dy);
+
+ // normal pressure contribution P_rhorho
+ for (ibin = 0; ibin < nbins; ibin++) {
+ // completely inside of R
+ if (rjsq < R2[ibin]) continue;
+
+ C=risq-R2[ibin];
+ D=B*B-4.0*A*C;
+
+ // completely outside of R
+ if (D < 0.0) continue;
+
+ sqrtD=sqrt(D);
+ alpha1=0.5*(-B+sqrtD)/A;
+ alpha2=0.5*(-B-sqrtD)/A;
+
+ if ((alpha1 > 0.0) && (alpha1 < 1.0)) {
+ zR=zi+alpha1*dz;
+ if ((zR < zhi) && (zR > zlo))
+ {
+ xR=xi+alpha1*dx;
+ yR=yi+alpha1*dy;
+ fn=fpair*fabs(xR*dx+yR*dy);
+
+ Pr_temp[ibin]+=fn;
+ }
+ }
+ if ((alpha2 > 0.0) && (alpha2 < 1.0)) {
+ zR=zi+alpha2*dz;
+ if ((zR < zhi) && (zR > zlo)) {
+ xR=xi+alpha2*dx;
+ yR=yi+alpha2*dy;
+ fn=fpair*fabs(xR*dx+yR*dy);
+
+ Pr_temp[ibin]+=fn;
+ }
+ }
+ }
+
+ // azimuthal pressure contribution (P_phiphi)
+ for (int iphi = 0; iphi < nphi; iphi++) {
+ alpha=(yi-xi*tangent[iphi])/(dx*tangent[iphi]-dy);
+
+ // no intersection with phi surface
+ if ((alpha >= 1.0) || (alpha <= 0.0)) continue;
+
+ // no contribution (outside of averaging region)
+ zL=zi+alpha*dz;
+ if ((zL > zhi) || (zL < zlo)) continue;
+
+ xL=xi+alpha*dx;
+ yL=yi+alpha*dy;
+
+ L2=xL*xL+yL*yL;
+
+ // no intersection (outside of Rmax)
+ if (L2 > R2kin[nbins-1]) continue;
+
+ ftphi=fabs(dx*ephi_x[iphi]+dy*ephi_y[iphi])*fpair;
+
+ // add to appropriate bin
+ for (ibin = 0; ibin < nbins; ibin++) if (L2 < R2kin[ibin]) {
+ Pphi_temp[ibin]+=ftphi;
+ break;
+ }
+ }
+
+ // z pressure contribution (P_zz)
+ for (int zbin = 0; zbin < nzbins; zbin++) {
+ // check if interaction contributes
+ if ((x[i][2] > binz[zbin]) && (x[j][2] > binz[zbin])) continue;
+ if ((x[i][2] < binz[zbin]) && (x[j][2] < binz[zbin])) continue;
+
+ alpha=(binz[zbin]-zi)/dz;
+
+ xL=xi+alpha*dx;
+ yL=yi+alpha*dy;
+
+ L2=xL*xL+yL*yL;
+
+ if (L2 > R2kin[nbins-1]) continue;
+
+ ftz=fabs(dz)*fpair;
+
+ // add to appropriate bin
+ for (ibin = 0; ibin < nbins; ibin++) if (L2 < R2kin[ibin]) {
+ Pz_temp[ibin]+=ftz;
+ break;
+ }
+ }
+ }
+ }
+
+ // calculate pressure (force over area)
+ for (ibin = 0; ibin < nbins; ibin++) {
+ Pr_temp[ibin]*=PrAinv[ibin]*Rinv[ibin];
+ Pphi_temp[ibin]*=PphiAinv;
+ Pz_temp[ibin]*=PzAinv[ibin];
+ }
+
+ // communicate these values across processors
+ MPI_Allreduce(Pr_temp,Pr_all,nbins,MPI_DOUBLE,MPI_SUM,world);
+ MPI_Allreduce(Pphi_temp,Pphi_all,nbins,MPI_DOUBLE,MPI_SUM,world);
+ MPI_Allreduce(Pz_temp,Pz_all,nbins,MPI_DOUBLE,MPI_SUM,world);
+
+ // populate array
+ for (ibin = 0; ibin < nbins; ibin++) {
+ array[ibin][0]=R[ibin];
+ array[ibin][2]=Pr_all[ibin]*nktv2p;
+ array[ibin][3]=Pphi_all[ibin]*nktv2p;
+ array[ibin][4]=Pz_all[ibin]*nktv2p;
+ }
+
+}
+
+/* ----------------------------------------------------------------------
+memory usage of data
+------------------------------------------------------------------------- */
+
+double ComputePressureCyl::memory_usage()
+{
+ double bytes =
+ (3.0*(double)nphi + 16.0*(double)nbins+5.0*(double)nbins) * sizeof(double);
+ return bytes;
+}
diff --git a/src/USER-MISC/compute_pressure_cylinder.h b/src/USER-MISC/compute_pressure_cylinder.h
new file mode 100644
index 0000000000000000000000000000000000000000..7151cca6cedab19f2439a139eb895b948f0aba86
--- /dev/null
+++ b/src/USER-MISC/compute_pressure_cylinder.h
@@ -0,0 +1,74 @@
+/* -*- 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.
+------------------------------------------------------------------------- */
+
+#ifdef COMPUTE_CLASS
+
+ComputeStyle(pressure/cylinder,ComputePressureCyl)
+
+#else
+
+#ifndef LMP_COMPUTE_PRESSURE_CYLINDER
+#define LMP_COMPUTE_PRESSURE_CYLINDER
+
+#include "compute.h"
+
+namespace LAMMPS_NS {
+
+class ComputePressureCyl : public Compute {
+ public:
+ ComputePressureCyl(class LAMMPS *, int, char **);
+ ~ComputePressureCyl();
+ void init();
+ void init_list(int, class NeighList *);
+ void compute_array();
+ double memory_usage();
+
+ private:
+ int nbins,nphi,nzbins;
+ double *Pr_temp,*Pr_all,*Pz_temp,*Pz_all,*Pphi_temp,*Pphi_all;
+ double *R,*Rinv,*R2,*PrAinv,*PzAinv,PphiAinv;
+ double Rmax,bin_width,nktv2p;
+ double *R2kin,*density_temp,*invVbin,*density_all;
+ double *tangent,*ephi_x,*ephi_y;
+ double *binz;
+
+ double zlo,zhi;
+
+ class NeighList *list;
+
+};
+
+}
+
+#endif
+#endif
+
+/* ERROR/WARNING messages:
+
+E: Illegal ... command
+
+Self-explanatory. Check the input script syntax and compare to the
+documentation for the command. You can use -echo screen as a
+command-line option when running LAMMPS to see the offending line.
+
+E: No pair style is defined for compute pressure/cylinder
+
+Self-explanatory.
+
+E: Pair style does not support compute pressure/cylinder
+
+The pair style does not have a single() function, so it can
+not be invoked by compute pressure/cylinder.
+
+*/
+