If you copy and paste the program source code below into a text editor, you can easily compile and run this.
(This class file name is hyios, so save this text editor as "hyios.java", and compile it.)
In this program, we first input the distance between two H atoms (= angstrom ) and exponential parameter Z, and the program outputs H2+ bond energy and antibond energy change from the original H atom and H+ proton separated by infinite distance.
Next input x-coordinate of the electron's position, and this program outputs the kinetic, Coulomb and total energies at this x position.
import java.util.Scanner;
class hyios {
public static void main(String[] args) {
// H2+ molecule, various Coulomb energy
Scanner stdIn=new Scanner(System.in);
System.out.println("distance between two H (= x angstrom )？ ");
double r=stdIn.nextDouble();
// input internuclear distance (= angstrom) of H2+
System.out.println("exponential parameter Z ？ ");
double z=stdIn.nextDouble();
r=r*z/0.5291772;
// r x Bohr radius = actual distance
// replacing Bohr radius by (bohr radius)/z
double me=9.1093826e-31; // me = electron's mass
double pai=3.141592653589793;
double epsi=8.85418781787346e-12; // epsi = permittivity
double h=6.62606896e-34; // h = Planck constant
double ele=1.60217653e-19; // ele = electron charge e
double aa=5.291772e-11/z; // aa = (bohr radius)/z
double ab=5.291772e-11;
double ea = Math.exp(-r); double eb = Math.exp(-2.0*r);
double ec=Math.exp(r);
double ss = ea * ( 1.0 + r + r*r/3.0 );
// ss = overlap integral
double sc = (1+r)*ea/aa;
// sc = K = overlap potential integral
// between nucleus-electron
sc = (-1.0 + (z-1.0))*(ele*ele)/(4.0*pai*epsi) * sc;
// (z-1.0) = exponent z (= besides 1 ) correction
double hc= -(ele*ele)/(4.0*pai*epsi) * ( 1.0/(r*aa) - eb*( 1.0/aa+ 1.0/(r*aa)) );
// hc = Coulomb attractive potential energy
// between elecron H-A wavefunction and the other nucleus-B
double hcc=((z-1.0)*ele*ele)/(aa*4.0*pai*epsi);
// hcc = exponent z (= besides 1 ) potential correction
double nucl=(ele*ele)/(4.0*pai*epsi*r*aa)*6.241509e18;
// nucl = Coulomb repulsive energy between two nuclei (eV)
hc=hc+hcc; hc=hc* 6.241509e18; // convert J to eV
sc=sc* 6.241509e18;
double total=-13.606*z*z+nucl+(hc+sc)/(1.0+ss);
double anti=-13.606*z*z+nucl+(hc-sc)/(1.0-ss);
double bondenergy = total-(-13.606);
double antibond = anti-(-13.606);
System.out.printf("H2+ bondenergy(eV):%.6f\n", bondenergy);
System.out.printf("H2+ antibond energy (eV):%.6f\n", antibond);
// calculation kinetic and potential energies at point x
r=r/z;
System.out.println(" ");
System.out.println("x-coordinate (= angstrom )？ ");
double x=stdIn.nextDouble();
x=x/0.5291772;
double ra= Math.sqrt((r*0.5+x)*(r*0.5+x));
// ra=distance between point-x and HA nucleus
double rb= Math.sqrt((r*0.5-x)*(r*0.5-x));
// rb=distance between point-x and HB nucleus
double orikine=13.606*2/rb-13.606;
// orikine = original kinetic energy (eV) at point x
// of the original H atom-B ( Z = 1 )
System.out.printf("original kinetic energy (eV):%.6f\n", orikine);
double era= Math.exp(-z*ra); double erb= Math.exp(-z*rb);
double nume= (-z*z+(2.0*z)/ra)*era+(-z*z+(2.0*z)/rb)*erb;
double kinex=(nume)/(era+erb)*13.606;
System.out.printf("kinetic energy-symme (eV):%.6f\n", kinex);
nume= (-z*z+(2.0*z)/ra)*era-(-z*z+(2.0*z)/rb)*erb;
double antikinex =(nume)/(era-erb)*13.606;
System.out.printf("kinetic energy-anti (eV):%.6f\n", antikinex);
double coulomb= -13.606*2/rb;
System.out.printf("original Coulomb energy (eV):%.6f\n", coulomb);
double tcoulo= coulomb-13.606*2*(1.0/ra-1.0/r);
System.out.printf("total Coulomb potential energy-x (eV):%.6f\n", tcoulo);
double symtotal= tcoulo + kinex + 13.606;
double antitotal= tcoulo + antikinex + 13.606;
// total energy change from the original isolated hydrogen
System.out.printf("total energy-symme (eV):%.6f\n", symtotal);
System.out.printf("total energy-antisym (eV):%.6f\n", antitotal);
}}