Hello,

We recently did some simulations to calculate the viscosity using the Muller-Plathe method. Using units=real, when we divide the momentum flux (j) over the velocity gradient, I thought to obtain the viscosity directly in Poise units (without any further conversion). However, this way I get crazy values for the viscosity.

Is there some conversion factor that I am missing?

Kind regards

Theo

Hello,

I guess the output of the fix, namely the total momentum, is in units of

(grams/mole) * (A/fs). Then one has to divide it by time (in fs) and

cross-sectional area (in A^2) to get the momentum flux. Finally the latter has

to be divided by the shear rate (in fs^-1). The result should therefore be in

units of (grams/mole)/(A*fs). To get the viscosity in Pa s = kg/(m s), one has

to apply the following conversion factor:

1 gram/mole/A/fs = 10^-3 / 6.02e23 kg * 10^10 m^-1 * 10^15 s^-1

I tried to use the Muller-Plathe method recently and to get results consistent

with other methods, I had to divide the results by 2, namely really apply the

definition (eq 3) of Muller-Plathe article. Is it right? In this case, maybe

this should be precised in the documentation...

Best regards,

Laurent Joly

Laurent is correct about the units - you need to convert mass*velocity

(grams/mole) * (A/fs) to Poise.

As the doc page for fix viscosity states:

As described below, the total momentum transferred by these velocity swaps is computed by the fix and >can be output. Dividing this quantity by time and the cross-sectional area of the simulation box yields a >momentum flux. The ratio of momentum flux to the slope of the shear velocity profile is the viscosity of >the fluid, in appopriate units. See the Muller-Plathe paper for details.

The cummulative momentum transferred between the bottom and middle of the simulation box (in the >pdim direction) is stored as a scalar quantity by this fix. This quantity is zeroed when the fix is defined and >accumlates thereafter, once every N steps. The units of the quantity are momentum = mass*velocity.

You also need to include the factor of 2 that the MP paper includes in

its formula, since the

system is relaxing in 2 directions due to periodicity. (I'll note

this in the doc page).

LAMMPS does not do these things for you, b/c you might use this fix in

different ways.

You might apply it to a subset of atoms, or use the flux quantity for

some other calculation.

If you want a viscosity, you need to use the correct formulas and

do the correct unit conversion.

Steve