I am interested in monitoring the electric field on a given atom in a charged system. I initially wrote some code to do this by computing the electric field from forcecoul in the compute() of function pair_lj_cut_coul_long.cpp. After seeing strange fluctuations in the magnitude of the electric field, and some further testing / debugging, I realized that compute() is not the right place to get this information. Basically the forces in compute() are not the same ones output using dump custom, so there are some steps between the Verlet iteration steps and end_of_step that do things I do not know.
I would like to get some advice from people who know the code better than I do on how to do this. I also need to know when (final_integrate, end_of_step, etc) is the best time to implement this. I do not want to redo the force calculation, so ideally I should be able to use the force calculation that is already done at every timestep. Basically I think I need to do something inside compute() since that's where forcecoul is, but I am missing what I need to do after that to get the correct (continuous) force and related electric field.
I would appreciate any advice on this matter. Thanks for your help!