(* K atom Kimball, Ne-center, outer-sp-spheres 21.12.2011/05.07.12 *)
Clear[k1,k2,k3,k4,sig1,sig2,sig3,sig4,c,z,R1,R2,R3,R4,S2,S3];
c = {k1 -> 1.0, k2 -> 1.0, k3 -> 1.0, k4 -> 1.0, sig1 -> 0.3, sig2 -> 0.3, sig3 -> 0.3, sig4 -> 0.3}; z=19.0;

T = 2.25*k1/R1^2+9.0*k2/R2^2 /. c;
ad = Sqrt[3./8.]; 
Vee=3.0*sig1/R1+12.0*sig2/R2+16/(R1+R2)+24*ad/(R1+R2) /. c;
Vne=-3.0*z/R1-8.0*z/(R1+R2);
S2 = R2*4^(1/3);

T = T + 9.0*k3/R3^2 /. c;
Vee = Vee+12.0*sig3/R3+80.0/(S2+R3)+24.0*ad/(S2+R3) /. c;
Vne = Vne-8.35*z/(S2+R3);
S3 = R3*2^(2/3);

T = T + 1.125*k4/R4^2 /. c;
Vee = Vee+18.0/(S3+R4);
Vne = Vne - z/(S3+R4);
S4=R4;

func = T + Vee + Vne;

t = FindMinimum[func, {R1,0.0721379}, {R2,0.244894}, {R3,0.7469272}, {R4,3.8494901},{Method -> Automatic},
{MaxIterations -> 500}]

N[-(Vee+Vne)/T /. c /. t[[2]],10]
S2 /. t[[2]]
S3 /. t[[2]]
S4 /. t[[2]]

plot1=Graphics3D[{Opacity[0.5],{Sphere[{0,0,0},S2],Sphere[{0,0,0},S3],
     Sphere[{0,0,0},S4]},{Darker[Green,1],Sphere[{0,0,0},R1]}}] /. t[[2]];
     
Show[plot1,{AspectRatio -> Automatic,Axes -> True}]