Hi,

I am wondering what is the unit for the stress tensor that is dumped by LAMMPS when using real units. Dump does not have units option as far as I know. units real is my command at the beginning of the input.

When I look at the website:

S = m.vi.vj + 0.5*SUM{ (ai-bi)Fij }

m in real units is given as g/cm3 in unit command, (ai-bi) is in angstroms and force is defined in kcal/mol.Angstroms in units command. v is in angstroms/fs.

When we divide by the volume later, volume is in angstroms cube.

These units are really not consistent, LAMMPS should be doing some conversions in the final output. I have found efactor, pfactor and boltz in initialize.h. But I am not sure in which units it writes in dump file for the stress tensor.

I have the similar unit problem with the velocities that are dumped. They really should be in angstroms/fs. But when I calculate temperature elsewhere using the initial velocities created by create gauss command to give 300.15K, I don’t get 300.15K but much high T.

When I calculate temperature, I use the following formula which is the same principle as for LAMMPS.

Ekin = ekfactor*(0.5*m*v^2)

T = 2*Ekin/(3*boltz)/n

I use eV for kinetic energy, A/fs for volume and amu for atom mass. K for temperature. And I have ekfactor for conversion which I calculated to be 103.64269. I calculate a much higher temperature for the same velocities but LAMMPS print 300K in the output.

I would appreciate if someone can clarify me the units in stress tensor and the velocity.

Thank you,

Burcu

Burcu,

From the dump command documentation here:

http://www.cs.sandia.gov/~sjplimp/lammps/doc/dump.html

We have:

"Note that this formula for stress does not include

virial contributions from intra-molecular interactions

(e.g. bonds, angles, torsions, etc). Also note that

this quantity is the negative of the per-atom pressure

tensor. It is also really a stress-volume formulation.

It would need to be divided by a per-atom volume to

have units of stress, but an individual atom's volume

is not easy to compute in a deformed solid.

Computation of stress tensor components requires a

loop thru the neighbor list and inter-processor

communication, so it can be inefficient to dump this

quantity too frequently or to have multiple dump

commands, each with stress tensor attributes."

Looking at the source code of fix_stress.cpp, we see:

// convert to pressure units (actually stress/volume =

pressure)

double nktv2p = force->nktv2p;

for (i = 0; i < nlocal; i++) {

stress[i][0] *= nktv2p;

Which means that LAMMPS is giving you units of

pressure*volume. For real units, this would be

atm*angs^3.

Mass has units of grams/mole (if units are real), not

gm/cm^3. See:

http://www.cs.sandia.gov/~sjplimp/lammps/doc/units.html

Note that nktv2p above does the units conversion.

My guess on the issue you're seeing with the

temperature and velocities not agreeing is that you

don't have the right value for n. If you're using

shake, degrees of freedom are removed and n is no

longer the number of atoms. Otherwise, the issue is

likely with the units conversion. If you're using real

units, volume would be angs^3, ke in kcal/mol, and T

in K. You'd have to use the proper conversion factors

like LAMMPS uses. For real units, you should have:

force->boltz = 0.001987191;

force->mvv2e = 48.88821 * 48.88821;

where the "boltz" is the Boltzmann factor, and mvv2e

converts from mass*velocity^2 to energy (in kcal/mol).

The conversion factors are in update.cpp.

Paul

I would appreciate if someone can clarify me the units in stress

tensor >and the velocity.

The doc page for the dump command doesn't state this explicity (I'll add it),

but all dump outputs are in the units described in the units doc pate. I.e.

if you dump velocities then they are in the same units as all velocities

that are input/output by LAMMPS, namely Ang/fmsec for real units or

Ang/psec for metal units, etc.

For the stress attributes (sxx, syy, etc), as the dump doc says, these

quantities are stress*volume, which means pressure * distance^3. So

for real units, it should be atmospheres * Ang^3.

Hope that helps,

Steve