(* H2+ sp-Molecule: H2+ two 0.5*(1p) overlapping 07.09.2012 *)
Clear[Q,R,d,c2,z,x];
c1=0.0; c2=(1-c1)/2; c22=2*c2; (* electron partition *)
x=-0.6;                        (* artificial separation, collapses when in min.fct. *)
Ekin=1.018*c22*9/(8*Q^2);
d=Q+R+x;
z=1.092;                       (* artificial eff. proton charge *)
Vnn=z^2/(2*d);
Vne=-c2*z*(3.0-(R/Q)^2)/Q-c22*z/(R+2*(Q+x));
If[Abs[R]>Q,Vne=-c22*z/(Q+Abs[R])-c22*z/(R+2*(Q+x))];
func=Ekin+Vnn+Vne;
t=FindMinimum[func,{Q,1.6},{R,0.349}]
u=t[[2]];
vir=-(Vnn+Vne)/Ekin /. u
hh=2*d /. u
force=z^2/hh^2-c2*z*(R/Q)^2/R^2-c2*z/(R+2*(Q+x))^2 /. u (* Hellmann-Feynman force not zero *)

plot1=Graphics[{
Circle[{-Q-x,0},Q],
Circle[{Q+x,0},Q],{Dashing[{0.02,0.02}],Circle[{Q+x,0},Abs[R]],Circle[{-Q-x,0},Abs[R]]},
Disk[{d,0},0.05],Disk[{-d,0},0.05],{Thickness[0.01],Line[{{-d,0},{d,0}}]}}] /. u;

Show[plot1,{AspectRatio->Automatic, Frame -> True, Axes -> True, GridLines -> Automatic,
PlotRange -> {{-2,2},{-1.5,1.5}}}]