Potentials from nonadiabatic dfn may not be correct

[David Dickinson]

Probably need to change from

phi = antot / sum(spec%dens*spec%z*spec%z / spec%temp)

to something like

‌

   g_work = 1.0
   phi = 0.0
   call integrate_species(g_work, spec%dens*spec%z*spec%z / spec%temp, phi)
   phi = antot / real(phi)

Worth double checking the other expressions in this routine.

[David Dickinson]

Also need to add in adiabatic contribution when needed.

[David Dickinson]

Adiabatic contribution added in PR #616. That also adds some test cases. Interestingly those tests pass even though they don’t implement the integrate_species approach. It turns out that when velocity resolution is high enough it seems that call integrate_species(g_work, spec%dens*spec%z*spec%z / spec%temp, phi) should tend to sum(spec%dens*spec%z*spec%z / spec%temp). This raises a question as to if we want to take advantage of the simplified equations for the nonadiabatic field solve to evaluate the integral analytically or if we want to stick with numerical integration to remain more consistency with the numerical errors elsewhere.

‌