Dear all,

I need to compute the average of the SQUARED energy

(Kinetic+Potential) over a certain interval of time in a dynamic

simulations:

<E_total^2>

I have seen how to compute the average of the kinetic and potential

energy independently.

But I couldn't find on the manual how to define variables on their sum

and squares.

Is it possible to define a quantity

E_tot = kinetic+potential

and then

E_square = E_tot*E_tot

and then use ave and reduce

?

What precise voice of the manual should I look at ?

How is the syntax ?

Any suggestions ?

In addition, how can I print the average of the distance between all

couples of atoms <r_ij> ?

And its square <r_ij^2> ?

Is it possible with ave and reduce ?

I have found just the command for the displacement from its original coordinates

compute displace/atom, but that is not what I need...

Thank you in advance for your help.

Bests,

Daniele Scopece

The first part is easy: You can define variables that do this, e.g.:

variable U equal pe

variable K equal ke

variable E equal “v_U + v_K”

variable E2 equal “v_E*v_E”

etc…

These can then be averaged with fix ave/time.

The second part is more difficult, and I am not sure if it is possible without modifications. There is a compute pair/local which gives distances between atoms involved in pairwise interactions, which you can dump or compute with dump/local and compute/local. Constructing the full NxN matrix of distances is more difficult though, and to me it seems that this is way easier in post-processing.

Compute reduce can take the output of compute pair/local (for distance)

and average it to a single value. You could also use fix ave/histo with

output from pair/local to calc a histogram of distances.

Neither can do those operations on the square of distances, and

there is no variable option that works for “local” values. I guess

either compute pair/local or compute reduce would need a new option

to create a squared quantity.

Steve

Keep in mind your r and r^2 descriptors for point 2 may turn out to scale with system size. At least that is the case of a homogeneous system (constant number density).

Radial intergration gives:

= \int r * P® dr = \int r * 4\pi*r^2/V_tot dr = 4\pi/V_tot * \int r^3 dr

using V_tot ~ r_max^3

one gets for ~ r_max a linear dependency with r_max being the maximum integration radius.

Carlos

PS: you may want to check my “fast food” approach to the calculation

Added a sumsq and avesq option to compute reduce.

Will be in the next patch.

Steve