Writing time-averaged per-chunk data (spread per-atom) to a dump-file

Dear LAMMPS experts,

I have a problem with sending to a dump the time-averaged per-chunk data.
Briefly, it can be described as follows:

  1. I compute a value per chunk that calculates a global array:
    compute stemp all temp/chunk molch temp adof 2

  2. Then I average over time the first (and the only) vector from the resulting array.
    It is noteworthy that the fix produces correct data in the dedicated file (tmp.out):
    fix aved1 all ave/time 1 1 1 c_stemp[1] file tmp.out mode vector

  3. Then I spread the result to a per-atom vector and write it to a dump-file:
    compute ctemp all chunk/spread/atom molch f_aved1
    dump 1 all custom 1 dump.out id type c_ctemp

For some reason, the answer is written only for the first chunk. For the rest chunks, only zeros are written in the dump-file.

I did my best to find the error and tried many ways to use wildcard and other LAMMPS scripting. I attach a minimal working example (input script and the input dump, which is a copy of lammps/examples/rdf-adf/data.spce) and I would be grateful for the help.

I expect it to be a realistic way to use LAMMPS since time-averaging of per-chunk data is mentioned in docs of chunk/spread/atom:
“… For inputs that are fixes, they should be a fix that calculates per-chunk values. For example, fix ave/chunk or fix ave/time (assuming it is time-averaging per-chunk data).”

Also, I would mention that another sequence (compute a value per chunk -> spread it over atoms -> average per-atom values over time -> dump) writes to dump NaNs and denormalized floats, which is probably a related issue.

Best regards,
Vsevolod Nikolskiy

units real
boundary p p p
atom_style full
bond_style harmonic
angle_style harmonic
atom_modify map array
pair_style lj/cut/tip4p/long 1 2 1 1 0.1546 9.0 9.0
suffix off
newton on
kspace_style pppm/tip4p 1.0e-5
suffix on

read_data data.spce

pair_coeff * * 0.0 0.0
pair_coeff 1 1 0.1852 3.1589
bond_coeff 1 0.0 0.9572
angle_coeff 1 0.0 104.52
group water type 1 2
fix 1 water shake 1.0e-4 200 0 b 1 a 1
fix 2 water nve
thermo 100
thermo_style custom step etotal ke pe temp evdwl ecoul elong press
velocity water create 300 123

#!!!
variable mol atom floor((id-1)/3+1)
compute molch all chunk/atom v_mol
compute stemp all temp/chunk molch temp adof 2
fix aved1 all ave/time 1 1 1 c_stemp[1] file tmp.out mode vector
compute ctemp all chunk/spread/atom molch f_aved1
dump 1 all custom 1 dump.out id type c_ctemp
#!!!
dump_modify 1 sort id
timestep 1
run 2

in.tip4p (1.19 KB)

data.spce (414 KB)

you are overlooking a subtle piece of information in the documentation of compute temp/chunk.
it computes 3 types of output:

  • a global scalar
  • a global vector of length 6
  • a global array of length nchunk with 1-3 columns

you are accessing the vector and hence you get the zeros or garbage data. the following corrects this and thus works as expected.

variable mol atom floor((id-1)/3+1)
compute molch all chunk/atom v_mol
compute stemp all temp/chunk molch temp adof 2
compute ctemp all chunk/spread/atom molch c_stemp[1]
fix ave all ave/atom 1 1 1 c_ctemp
dump 1 all custom 1 dump.lammpstrj id type v_mol f_ave c_ctemp
dump_modify 1 sort id

timestep 1
run 2