Local Stress Tensor

[RudolfWeeber]

On the mailing list, interest was expressed in the Irving Kirkwood local stress tensor: https://espressomd.github.io/doc4.0.2/analysis.html#local-stress-tensor The implmentation was rmoved in 4.1, since it was incorrect.

A reimplementation using current infrastructure might be around the following lines

• A histogram with Vector3i as index type and Vector9d as value type
  • The histogram needs a ghost layer, i.e., not only the volume of the MPI rank has to be covered, but also one interaction range beyond that. The reason is that each particle pair is visited only once. In particular, if two particles reside on different mpi ranks, only one of the ranks calculates the force, which is added to both particles (one of which is a ghost for the particle on the other rank). Th force on the ghost is transferred onto its host particle later during ghost force reduciton.
  • A parallel reduciton for the ghost layers of the histogram is not needed: the histograms of each mpi rank need to be sent to the head node anyway, to return them to the Python interface. Hence, the head node can add up values from the ghost layers
• Use the short range loop (which executes a kernel on non-bonded interacting particle pairs and on pairs of particles connected by a bond
  • Calcualte the pair force (as in the normal force calculation)
  • Calculate the pressure contribution of the pair using outer_product(distance,force) (as in normal virial pressure calculaiton)
  • Wirte std::pair<Vector3i, double> get_line_bin_intersection() which calculates bin indeces and wieghts for all histogram binds intersected by the line connecting the particles. Make sure to handle periodic boudnary conditons correctly: line from p1.pos() to p1.pos() +minimal_image_distance(p2,p1,box_geo) or the sth of that kind
  • Use a loop over the result of that function to add parts of the pair pressure contribution to the resepctive histogram bins

Sample code for using the short_range_loop can be foudn in forces.cpp:force_calc() The pair virial code is in pressure_inline.hpp: add_nonbonded_pair_virials() and add_bonded_pressure_tensor()