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 remass, so save this text editor as "remass.java", and compile it.)
In this program, we first input the number ( 1-3 ) of reduced mass conditions.
In "Normal" condition (= 1 ), we use the reduced mass except when the center of mass of two electrons is just at the nucleus (= use reduced mass except when at the start where two electrons just start from x-axis which uses ordinary electron mass because the helium nucleus transiently stops only at the start ).
In "Not reduced mass" (= 2 ) condition, we always use the ordinary electron mass.
In "All reduced mass" (= 3 ) condition, we always use the reduced mass even when the helium nucleus transiently stops
In this program, first we input the initial x-coordinate r1 (in MM) of electron 1 (= the electron-1 and electron-2 start from x-axis ), and the absolute value of the total energy E (in eV) of helium atom.
↑ This initial input r1 value increases per calculation of the electron's 1/4 orbit until r1 + 30.
From the inputted values, this program computes electrons' 1/4 orbits (= electron-1 moves from x-axis to y-axis on x-y place, while electron-2 moves from x-axis to z-axis on the perpendicular x-z plane symmetrically ) by rigorously calculating Coulomb forces among two electrons and a nucleus, de Broglie wavelength, updating electrons' positions, velocities at short time intervals of 1SS (= 1 × 10-23 second ).
Then, this program outputs the x (= VX ) and y (= VY ) components of electron-1's last velocity after the electron has moved a quarter of its orbits from x-axis to y-axis (= electron-1 must cross the y-axis perpendicularly with the last VY velocity being zero to form the most stable symmetrical orbit around the nucleus = so we should choose the uniquely-determined r1 which gives the zero last VY in each different helium total energy E ).
This program also outputs WN (= the number of de Broglie's waves included in one quarter of the orbital, which should be 0.250000 when one orbit is 1 × de Broglie wavelength, which case uniquely gives the predicted helium energy by this real helium model's program ).
Here we use new convenient units of time = 1 SS = 1 × 10-23 second, distance = 1 MM = 1 × 10-14, velocity = 1MM/SS = 1 × 109 m/s, acceleration = 1MM/SS^2 = 1032 m/s^2.
The initial x-coordinate is automatically increased per calculation until +30 (= these different initial x-coordinates give different orbits, out of which we can uniquely determine only one initial x-coordinate that gives a symmetrical orbit around a nucleus with zero y-component of electron-1's last velocity (= VY = 0.000000 ) after electron-1 moves a quarter of its orbit ).
import java.util.Scanner;
class remass {
public static void main(String[] args) {
// use reduced mass or electron mass ?
// electron-1 = e1
// input initial e1's x-coordinate = r1 (MM)
Scanner stdIn=new Scanner(System.in);
/*
choose "always use reduce mass = 3",
"always use ordinary electron mass = 2"
"Normal = 1 = use reduced mass except when two electrons are just at the opposite sides of the nucleus (= use ordinary elecrton mass only at starting point where electrons-1 and 2 just start from x-axis )"
*/
System.out.println("Input 1 (=Normal), 2 (=NOT reduced mass) or 3 (=All reduced mass)");
double ree=stdIn.nextDouble();
System.out.println("r1 between nucleus and electron 1 (MM)？ ");
double r=stdIn.nextDouble();
// input absolute value of helium energy = |E| eV
System.out.println("absolute value of total energy |E| of helium atom (eV) ? ");
double E=stdIn.nextDouble();
// me = electron mass kg
double me=9.1093826e-31;
// epsi = electric constant, pai = pi
double pai=3.141592653589793;
double epsi=8.85418781787346e-12;
// h = Planck constant, ele = electron charge
double h=6.62606896e-34;
double ele=1.60217653e-19;
// Z = helium atomic number
double Z = 2.0;
// nucle = helium nuclear mass kg
double nucle = 6.64465650e-27;
// rm = reduced mass
double rm=(2.0*me*nucle)/(2.0*me+nucle);
rm=rm*0.5;
// when Not using reduced mass
if ( ree > 1.5 && ree < 2.5 ) { rm = me ;}
// repeat until r1 increases to r1+30
for (int i=1;i < 30;i++) {
// poten = initial potential energy (J)
double poten=-(2.0*Z*ele*ele)/(4.0*pai*epsi*r*1.0e-14)+(ele*ele)/(4.0*pai*epsi*2.0*r*1.0e-14);
//vya= total E ( eV to J )-potential energy
double vya=-(E*1.60217646e-19)-poten;
if (vya > 0) {
// vyb=electron initial velocity (m/sec)
double vyb=Math.sqrt(vya/me);
// when always using reduced mass
if ( ree > 2.5 ) { vyb=Math.sqrt(vya/rm); }
// electron-1 = e1
// change m/sec (= vyb ) to MM/SS (= VY)
// 1 MM/SS = 10^9 m/s
// VX, VY is initial velocity of electron-1
double VY=vyb*1.0e-9; double VX = 0.0;
// prexx, preyy = initial e1's x,y-coordinate
double prexx=r; double preyy=0.0;
// WN = number of de Broglie wavelength
double WN=0.0;
double yy,vk, leng,wav; double xx=0.0;
do {
// xx, yy = electron 1 position after 1SS
xx=prexx+VX; yy=preyy+VY;
vk=VX*VX+VY*VY;
// leng = electron's moving length (m) for 1SS
leng=Math.sqrt(vk)*1.0e-14;
// wav = de Broglie wavelength = h/mv (m)
wav=h/(rm*Math.sqrt(vk)*1.0e9);
// WN = sum of de Broglie wavelengths
// in all short segments
WN=WN+leng/wav;
// update e1's velocity from Coulomb force
// ra = distance between e1 and nucleus
double ra=Math.sqrt(prexx*prexx+preyy*preyy);
// rb = distance between two electrons-1, 2
double rb=Math.sqrt(4.0*prexx*prexx+2.0*preyy*preyy);
// change unit from MM into meter
ra=ra*1.0e-14; rb=rb*1.0e-14;
prexx=prexx*1.0e-14; preyy=preyy*1.0e-14;
double ac=(ele*ele)/(4.0*pai*epsi*rm);
// calculate accelertion (MM/SS^2 = 10^32 m/s^2)
// update e1's x-velocity VX from Coulomb
VX=VX+1.0e-32*ac*prexx*(-Z/(ra*ra*ra)+2.0/(rb*rb*rb));
// update e1's y-velocity VY from Coulomb
VY=VY+1.0e-32*ac*preyy*(-Z/(ra*ra*ra)+1.0/(rb*rb*rb));
prexx=xx;preyy=yy;
} while (xx >=0);
// repeat until e-1 reaches y-axis
// until electron-1 has moved 1/4 orbit
// display initial e1's coordinate = r1
System.out.print("r1: "+r+" ");
// display e1's last velocity = VX, VY
System.out.printf("VX:%.6f ", VX);
System.out.printf("VY:%.6f ", VY);
// WN = de Broglie wavelength of 1/4 orbit
System.out.printf("WN:%.6f\n", WN);
} r=r+1;
}}}