how to output an average of compute gyration/chunk

_Frischknecht_Amalie · September 1, 2017, 10:26pm

Hi,

I am running simulations of a simple polymer melt. I want to calculate the average radius of gyration of all the polymer chains at a given timestep. I tried to do this by computing the Rg for each molecule using gyration/chunk, and then averaging all of the resulting Rg’s to get one scalar variable. I think I have to use a variable to do this, since the output of compute gyration/chunk is a global vector. So I used the ave() function. But now I can’t seem to print the result to a file, using fix ave/time. Here is what I have tried:

compute c1 all chunk/atom molecule

compute myChunk all gyration/chunk c1

variable aveRg vector ave(myChunk)

fix 3 all ave/time 100 1 100 v_aveRg file Rg.out mode vector

But I get an error saying:

ERROR: Invalid special function in variable formula

If I instead tell fix ave/time that its mode is scalar, I get:

ERROR: Fix ave/time variable is not equal-style variable

But ave() requires a vector-style variable. So I’m stuck. Is this calculation possible in LAMMPS? Or should I just output the Rg of each individual molecule and use a post-processing step to average them?

Thanks,

Amalie

akohlmey · September 1, 2017, 10:44pm

Hi,

I am running simulations of a simple polymer melt. I want to calculate
the average radius of gyration of all the polymer chains at a given
timestep. I tried to do this by computing the Rg for each molecule using
gyration/chunk, and then averaging all of the resulting Rg’s to get one
scalar variable. I think I have to use a variable to do this, since the
output of compute gyration/chunk is a global vector. So I used the ave()
function. But now I can’t seem to print the result to a file, using fix
ave/time. Here is what I have tried:

compute c1 all chunk/atom molecule

compute myChunk all gyration/chunk c1

variable aveRg vector ave(myChunk)

shouldn't aveRg be an equal style variable? the output of ave() is a
scalar, right?

axel.

sjplimp · September 6, 2017, 2:42pm

The confusion may be that this message:

ERROR: Fix ave/time variable is not equal-style variable

refers to the aveRG variable, not the ave() function

or its arg.

Ave() takes a vector as an arg and produces a scalar.
So, as Axel says you need

to assign it to an equal-style variable, not vector-style.

Then you should be able to use mode scalar for fix ave/time.

Steve

_Frischknecht_Amalie · September 12, 2017, 2:53pm

Thanks, but I don’t think that solves my problem. If I instead set aveRg to be an equal style variable, like this:

compute c1 all chunk/atom molecule

compute myChunk all gyration/chunk c1

variable aveRg equal ave(c_myChunk)

fix 3 all ave/time 100 1 100 v_aveRg file Rg.out mode scalar

then I get the following run-time error:

ERROR: Invalid special function in variable formula (…/variable.cpp:4015)

My understanding is that the gyration/chunk compute returns a global vector, as per the lammps manual: “This compute calculates a global vector if the tensor keyword is not specified and a global array if it is.”

So the input to ave() is a vector, so I need to specify it with vector. The ouput of the ave() function should be a scalar. So I think I should use:

compute c1 all chunk/atom molecule

compute myChunk all gyration/chunk c1

variable aveRg vector ave(c_myChunk)

fix 3 all ave/time 100 1 100 v_aveRg file Rg.out mode scalar

But now LAMMPS complains with:

ERROR: Fix ave/time variable is not equal-style variable (…/fix_ave_time.cpp:216)

So I guess fix ave/time can’t deal with a function that takes a vector and returns a scalar? Any suggestions on how to access the result of the ave() function when it’s input is a vector?

Thanks,

Amalie

akohlmey · September 12, 2017, 3:02pm

Thanks, but I don’t think that solves my problem. If I instead set aveRg
to be an equal style variable, like this:

compute c1 all chunk/atom molecule

compute myChunk all gyration/chunk c1

variable aveRg equal ave(c_myChunk)

fix 3 all ave/time 100 1 100 v_aveRg file
Rg.out mode scalar

then I get the following run-time error:

ERROR: Invalid special function in variable formula (../variable.cpp:4015)

are you using an up-to-date version of LAMMPS?

if i add the lines below to the micelle input example right before the run
command, it works for me without a hitch:

compute c1 all chunk/atom molecule
compute myChunk all gyration/chunk c1
variable aveRg equal ave(c_myChunk)
fix 0 all print 1000 "aveRg = ${aveRg}"
fix xxx all ave/time 100 1 100 v_aveRg file Rg.out mode scalar

this is with LAMMPS (1 Sep 2017)

i get output like this:

Setting up Verlet run ...
  Unit style : lj
  Current step : 0
  Time step : 0.005
Per MPI rank memory allocation (min/avg/max) = 4.149 | 4.149 | 4.149 Mbytes
Step Temp E_pair E_mol TotEng Press
       0 0.45 -1.7056163 0.08808163 -1.1679097 3.9431686
aveRg = 1.0562786187384
    1000 0.45 -1.9727611 0.058607365 -1.4645287 1.9982388
aveRg = 1.0588691476084
    2000 0.45313774 -2.0077426 0.060108702 -1.4948738 1.5384533
aveRg = 1.06115831186531
    3000 0.45 -2.0244054 0.058635034 -1.5161454 1.5373956
aveRg = 1.05915140277307
    4000 0.4528704 -2.0263317 0.066702425 -1.5071363 1.3543361
aveRg = 1.073740217242
    5000 0.44407998 -2.0218661 0.063798656 -1.5143576 1.4110054
aveRg = 1.07419524452449
    6000 0.45385994 -2.0407844 0.055460161 -1.5318425 1.1077206
aveRg = 1.07896896871454
    7000 0.43736612 -2.0290443 0.061800741 -1.5302419 1.27996
aveRg = 1.07368344655849
    8000 0.44536091 -2.0617669 0.064805485 -1.5519716 1.103032
aveRg = 1.07357485478261
    9000 0.45321257 -2.0463507 0.0650004 -1.5285154 1.1931682
aveRg = 1.07122062028487
   10000 0.45828949 -2.0635306 0.061828977 -1.543794 1.0504713
aveRg = 1.06480897337791

My understanding is that the gyration/chunk compute returns a global
vector, as per the lammps manual: “This compute calculates a global vector
if the *tensor* keyword is not specified and a global array if it is.”

So the input to ave() is a vector, so I need to specify it with vector.
The ouput of the ave() function should be a scalar. So I think I should use:

compute c1 all chunk/atom molecule

compute myChunk all gyration/chunk c1

variable aveRg vector ave(c_myChunk)

fix 3 all ave/time 100 1 100 v_aveRg file
Rg.out mode scalar

But now LAMMPS complains with:

ERROR: Fix ave/time variable is not equal-style variable
(../fix_ave_time.cpp:216)

So I guess fix ave/time can’t deal with a function that takes a vector and
returns a scalar? Any suggestions on how to access the result of the ave()
function when it’s input is a vector?

my guess is that you have an old version of LAMMPS with a bug that has
since been fixed.

axel.

_Frischknecht_Amalie · September 12, 2017, 3:17pm

I was using the 11 Aug 2017 version of lammps. I tried the micelle example you list below and it worked fine. I then modified my script and it is now working as well—I think I had a couple of small errors before (at one point forgot the “c” to refer to “c_myChunk”) and didn’t try all combinations once I had fixed a couple of them. This gives me what I want:

compute c1 all chunk/atom molecule

compute myChunk all gyration/chunk c1

variable aveRg equal ave(c_myChunk)

fix 3 all ave/time 100 1 100 v_aveRg file Rg.out mode scalar

Thanks much!!

Amalie