[lammps-users] issue writing a compute temp biased with angular velocity

Dear all,

I am trying to write a compute temp/omega measuring the temperature of
a group after the center of mass velocity and the angular velocity
have been removed. The modified code seems to be working both in
serial and parallel when I only try to compute the temperature, and it
also works in serial when I try to bias my thermostat with compute
temp/omega. But in parallel the code hangs just after thermo output
for timestep 0 has been printed on the screen. I'm not at all an
expert in programming and here I don't have a clue on how to debug
this... Could someone help me on this?

Enclosed are the modified files with regard to current lammps
distribution, and the input file generating the problem...


torque.tgz (11.6 KB)

Dear all,

I finally found out that my issue came from using group functions (for
instance group->vcm) inside remove_bias and restore_bias, where the
MPI_Allreduce command hangs. This could be solved by calling
remove_bias_all and restore_bias_all in the thermostat.

However I noticed that remove_bias_all and restore_bias_all, though
present in each compute_*.cpp code, is not called by any thermostat in
the current lammps version. Only remove_bias and restore_bias are

So I had two questions:

* why did you chose at some point to only call remove_bias and
restore_bias, i.e. remove and restore bias atom by atom, instead of
calling remove_bias_all and restore_bias_all?

* If I slightly modify fix_nh.cpp to call remove_bias_all and
restore_bias_all (see below), would it work with all the other
computes (since the corresponding code exists in all others
compute_*.cpp)? If I understood correctly I believe it should provide
the exact same results, but I just wanted to check...


The bias_all() methods are called by other computes,
like compute temp/sphere. The thermostats use
the per-particle bias methods, so I wouldn't change
that. The idea is that the compute, say compute temp/rotate
in your case, is supposed to do any global operations it
needs to calculate the bias, when it's compute_scalar()
method is invoked to calculate the temperature. E.g. the
group operations you mention. This can involve creating
some data structure with the bias values if needed. Then
when remove/restore bias is called on a per-particle basis
it can do what is needed. All the thermostats call compute_scalar()
before they call remove_bias(), so this works. See
compute temp/com for an example.

Is there some reason that model doesn't work for what you
want to do?


Thanks a lot for this clarification! Doing as you wrote seems to work
like a charm. Just to check: in order to remove the angular velocity,
I have to store a velocity bias for all atoms, so I define: "double
**vbiasall;" and "int maxbias;" in "compute_temp_rotate.h".

Then I just have to add to ComputeTempRotate():

maxbias = 0;

and before I affect vbiasall[i][0,1,2] in compute_scalar() AND
compute_vector() (or only in compute_scalar?):

if (nlocal > maxbias) {
  maxbias = atom->nmax;
  vbiasall = memory->create_2d_double_array(maxbias,3,"temp/rotate:vbiasall");

(and destroy it at the end, and add the corresponding memory_usage() function)

Is that right?

Enclosed the whole source files.


compute_temp_rotate.cpp (8.31 KB)

compute_temp_rotate.h (1.35 KB)