Forces: Include ion-ion force contribution to ACForce

[rcclay]

Currently, ACForce does not appear to accumulate the ion-ion contribution to the total force on the ions. I feel like this should be included by default. Verification that this can be computed is needed for the PBC case, and implementation is definitely needed for the openBC case.

[jtkrogel]

@rcclay We are targeting molecules in open boundary conditions first, so ion-ion force in OBC will be needed. It seems like this should be pretty simple to add. Is this something you can do in the near future?

[rcclay]

@jtkrogel Sure. It's the easiest force to compute, just need to sit down, implement, and unit test it.

[tiihonej]

I can take a look at the ion-ion force in OBC, unless you already have it ready, @rcclay?

[rcclay]

Feel free to take a stab at it. It'll go in src/QMCHamiltonians/CoulombPotential.h . If you check the constructor, there's a "if(!is_active)" block, indicating whether QMCPACK expects those particles to be moving or not. If they don't, it currently computes the like-species energy in the constructor and returns it every time the estimator is queried so you don't have to recompute it at every electron step. I would do forces the same way. You can consult src/QMCHamiltonians/CoulombPBCAA.cpp to see how this was done in PBC's. OpenBC's will be easier since you don't have to worry about ewald breakups.

1.) Implement bare like species forces. Look at evaluateAA() in CoulombPotential to get an idea of how the charge and distance tables are accessed. You'll be messing with \mathbf{dr} now instead of just scalar quantities for forces. You can park it in its own evaluateForces routine or something for compartmentalization.
2.) ForceBase class provides a prefix variable (the name of the force quantity in scalar.dat), an addObservablesF() function for registering these to scalar.dat, and a "forces" vector for storing the forces. You'll need to call the addObservablesF() and set the aforementioned variables, ideally with some logic guarding these operations (if ComputeForces==True, then do stuff... otherwise don't). See CoulombPBCAA
3.) To work with the ACForce estimator, you'll need to add an evaluateWithIonDerivs() function. You'll need to set "hf_terms" to whatever forces you computed, and set "pulay_terms" to zero. Likewise, see CoulombPBCAA.

Feel free to ask questions!

[rcclay]

Closed by #1820