Back to the chemical bond page
Top page ( correct Bohr model including the two electron atom )
This program is a little long. So if you copy and paste the below program source code into a text editor,
you can easily compile and run this.
This program's class file name is watr, so save this text editor as "watr.java", and compile it.
In this program, nuclei are gray circles.
Here we use the new units, ( 1 MM = 10-14 meter).
Each coordinate of electrons (+X (MM), +Y (MM), +Z (MM)) in the text box means "relative" position from these nuclei.
(ele 0-5 are from oxygen nucleus, ele 6 is from hydrogen nucleus 0 (H0), and ele7 is from H1.)
You can change the coordinate (+X, +Y, +Z) freely.
(Enter the values into the textboxes, and press the Enter key.)
"nuc (MM)" means the distance between these nuclei and electrons.
V (eV) and T (eV) means the potential, and kinetic energies of each electron.
tV (eV) is the total potential energy.
CF of ele0-5 means the force toward the center, CF of ele6,7 means the force toward H-O-H plane.
(fx, fy, fz) mean force components other than CF.
(FX, FY, FZ) mean the force component acting on each nucleus.
H0, On mean the force toward H0 nucleus and O nucleus.
Waves (wn) is number of de Broglie's waves contained in one orbit.
when you choose the O-H bond length (MM) in the scrollbar and click "O-H (MM)" button, the internuclear distance of O-H change.
And when you choose the angle in the scrollbar and click "angle" button, H-O-H angle changes.
import java.awt.*;
import java.awt.event.*;
import javax.swing.*;
import java.util.Scanner;
public class watr extends JPanel // virial theorem of H2O(water)
{
public static void main(String arg[])
{
JFrame frame = new JFrame("H2O (water)");
J2DPanel j2dpanel = new J2DPanel();
frame.getContentPane().add(j2dpanel); frame.setSize(1180,700);
frame.setVisible(true); frame.setBackground(Color.black);
frame.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
}
}
class J2DPanel extends JPanel implements ActionListener
{
double pai=3.141592653589793; double epsi=8.85418781787346e-12;
double h=6.62606896e-34; double elc=1.60217653e-19;
double me=9.1093826e-31; double suh=5200.0*5200.0; // suh=(Bohr radius)^2
JTextField elp[][]=new JTextField[8][11]; // elp=each electron's parameters
JTextField impho=new JTextField(7); // impho= total V (eV)
JTextField mmpho[][]=new JTextField[3][4]; // mmpho=each nuclear parameters
JButton b1=new JButton("O-H (MM)"); JButton b2=new JButton("angle");
String ope[]={"9000","9584","10000"}; // scrollbar of O-H distances
JComboBox coom=new JComboBox(ope);
String ope2[]={"90.0","104.45","110.0"}; // scrollbar of H-O-H angles
JComboBox coom2=new JComboBox(ope2);
double rtw=Math.sqrt(2); double rth=Math.sqrt(3);
double rsi=Math.sqrt(6); double rfi=Math.sqrt(5);
double hpr[][]=new double[8][11]; double hpr2[][]=new double[8][11];
double den=6.273; // den=central charge
double lengt=9584.0; double angl=104.45; // lengt=O-H distance, angl=H-O-H angle
double lengt2=lengt*Math.sin(angl*0.5*pai/180.0);
double lengt3=lengt*Math.cos(angl*0.5*pai/180.0);
//nux[n][0-2]=coordinate, n=1,2,3: O , H0 , and H1 nuclei
double nux[][]={{7857.00, 12500.00, 12500.00},
{7857.0+lengt3, 12500.00-lengt2, 12500.00},
{7857.0+lengt3, 12500.00+lengt2, 12500.00}};
// te1=initial conditions of each electron
double te1[][]={{3151,-3244, 0.0}, {3151, 3244, 0.0},
{-3290, 3290, 0}, {-3290, -3290, 0},
{0, 0, -4557}, {0,0, 4556},
{-260.0, 319.0, -4535.0}, {-260.0, -319.0, 4535.0}};
double te2[][]={{3214,-3215, 0}, {3214, 3215, 0.0},
{-3284, 3293, -10}, {-3284, -3293, 10},
{0, 0, -4562}, {0,0, 4562},
{-162.0, 330.0, -4435.0}, {-162.0, -330.0, 4435.0}};
double te3[][]={{3138,-3247, 0.0}, {3138, 3247, 0.0},
{-3290, 3290, -20}, {-3290, -3290, 20},
{0, 0, -4556}, {0,0, 4556},
{-260.0, 330.0, -4650.0}, {-260.0, -330.0, 4650.0}};
double te4[][]={{3120,-3240, 0.0}, {3120, 3240, 0.0},
{-3290, 3290, -20}, {-3290, -3290, 20},
{0, 0, -4557}, {0,0, 4557},
{-260.0, 310.0, -4310.0}, {-260.0, -310.0, 4310.0}};
double te5[][]={{3150,-3245, 0.0}, {3150, 3245, 0.0},
{-3290, 3290, -20}, {-3290, -3290, 20},
{0, 0, -4558}, {0,0, 4558},
{-200.0, 300.0, -4800.0}, {-200.0, -300.0, 4800.0}};
public J2DPanel()
{
setBackground(Color.black);
JPanel p=new JPanel();
p.setLayout(new GridLayout(11,12));
int aaa=0;
for (int el=0; el <=7; el++) { // elp[0-7][0-2]=textboxes of each electron's coordinate
for (int pos=0; pos <=2; pos++) {
elp[el][pos]=new JTextField(7); elp[el][pos].addActionListener(this);
hpr[el][pos]=0.0;
}}
for (int el=0; el <=7; el++) { // elp[0-7][3-10]=textboxes of other parameters
for (int pos=3; pos <=10; pos++) {
elp[el][pos]=new JTextField(7);
hpr[el][pos]=0.0;
}}
for (int el=0; el <=2; el++) {
for (int pos=0; pos <=3; pos++) {
mmpho[el][pos]=new JTextField(7);
}}
// layout
String sihy[]={"eNo ", "+X(MM)", "+Y(MM)", "+Z(MM)", "nuc(MM)",
"V(eV)", "T(eV)", "Force", "fx ", "fy", "fz", "Waves"};
for (int el=0; el <=11; el++) {
p.add(new Label(sihy[el]));
}
for (int el=0; el <=7; el++) {
p.add(new Label("ele "+el+" "));
for (int pos=0; pos <=10; pos++) {
p.add(elp[el][pos]);
}}
p.add(new Label("O nuc "));
for (int pos=0; pos <=3; pos++) {
p.add(mmpho[0][pos]);
} p.add(impho);
p.add(new Label("H0nuc "));
for (int pos=0; pos <=3; pos++) {
p.add(mmpho[1][pos]);
} p.add(new Label(" -- "));
p.add(new Label("H1nuc "));
for (int pos=0; pos <=3; pos++) {
p.add(mmpho[2][pos]);
} p.add(new Label(" -- "));
p.add(b1); p.add(coom); p.add(b2); p.add(coom2);
for (int pos=0; pos <=1; pos++) {
p.add(new Label(" -- "));
}
coom.setSelectedItem("9584"); b1.addActionListener(this);
coom2.setSelectedItem("104.45"); b2.addActionListener(this);
add(p,"South");
for (int el=0; el <=7; el++) {
double nnuc=Math.sqrt(te1[el][0]*te1[el][0]+te1[el][1]*te1[el][1]+te1[el][2]*te1[el][2]);
aaa=(int)(nnuc);
elp[el][3].setText(Integer.toString(aaa)); // show distance between nuclei and electrons
for (int jou=0; jou <=2; jou++) { // hpr[0-7][0-2]=each electron's coordinate
hpr[el][jou]=te1[el][jou];
if (el < 6) {hpr[el][jou]=hpr[el][jou]+nux[0][jou];}
if (el==6) {hpr[el][jou]=hpr[el][jou]+nux[1][jou];}
if (el==7) {hpr[el][jou]=hpr[el][jou]+nux[2][jou];}
elp[el][jou].setText(Integer.toString((int)(te1[el][jou])));
}}
} // public J2DPanel() end
public void actionPerformed(ActionEvent e) {
String ss;
int RR=0; double Rf1,Rf2,Rf3; int teis=0;
if (e.getSource() == b1) { RR=1; // when internuclear distance change (b1 click)
ss=(String)coom.getSelectedItem();
if (ss=="9000") {lengt=9000;} if (ss=="9584") {lengt=9584;}
if (ss=="10000") {lengt=10000;}
}
if (e.getSource() == b2) { RR=1; // when H-O-H angle change (b2 click)
ss=(String)coom2.getSelectedItem();
if (ss=="90.0") {angl=90.0;} if (ss=="104.45") {angl=104.45;}
if (ss=="110.0") {angl=110.0;}
}
if (RR==1) {
lengt2=lengt*Math.sin(angl*0.5*pai/180.0); lengt3=lengt*Math.cos(angl*0.5*pai/180.0);
double nuux[][]={{7857.00, 12500.00, 12500.00},
{7857.0+lengt3, 12500.00-lengt2, 12500.00},
{7857.0+lengt3, 12500.00+lengt2, 12500.00}};
for (int ett=0; ett <=2; ett++) {
for (int sws=0; sws <=2; sws++) {
nux[ett][sws]=nuux[ett][sws];
}}
for (int ett=0; ett <=7; ett++) { // each electron's coordinate reset
for (int sws=0; sws <=2; sws++) {
Rf1=te1[ett][sws];
if (lengt==9584 && angl==90.0) {Rf1=te2[ett][sws];}
if (lengt==9584 && angl==110.0) {Rf1=te3[ett][sws];}
if (lengt==9000 && angl==104.45) {Rf1=te4[ett][sws];}
if (lengt==10000 && angl==104.45) {Rf1=te5[ett][sws];}
elp[ett][sws].setText(Integer.toString((int)Rf1));
}}
}
for (int ett=0; ett <=7; ett++) { // when electron's positions change
for (int sws=0; sws <=2; sws++) {
ss=elp[ett][sws].getText(); Rf1=Double.parseDouble(ss);
Rf2=0.0;
// change relative coordinate to absolute coordinate
if (ett < 6) {Rf2=Rf1+nux[0][sws];}
if (ett==6) {Rf2=Rf1+nux[1][sws];}
if (ett==7) {Rf2=Rf1+nux[2][sws];}
Rf3=Rf2/71.428; // change MM to pixel
// upper and lower limits of coordinates
if (Rf3 > 349 && sws==0) {Rf3=349; Rf2=Rf3*71.428; teis=1;}
if (Rf3 > 349 && sws==2) {Rf3=349; Rf2=Rf3*71.428; teis=1;}
if (Rf3 > 349 && sws==1) {Rf3=349; Rf2=Rf3*71.428; teis=1;}
if (Rf3 < 1) {Rf3=1; Rf2=Rf3*71.428; teis=1;}
hpr[ett][sws]=Rf2;
if (ett < 6) {Rf1=Rf2-nux[0][sws];}
if (ett==6) {Rf1=Rf2-nux[1][sws];}
if (ett==7) {Rf1=Rf2-nux[2][sws];}
if (teis==1) {elp[ett][sws].setText(Integer.toString((int)(Rf1)));}
}}
repaint();
}
public void update(Graphics g)
{
paint(g);
}
public void paintComponent(Graphics g)
{
double kro,krr,krk,kro2,krr2,krk2,pot,pota,potb,pot2,pota2,potb2,
gx,gy,gz,ggx,ggy,ggz,ttav,toav;
kro=0.0; krr=0.0; krk=0.0; kro2=0.0; krr2=0.0; krk2=0.0;
int ex,ey,ez,xk,yk,zk; String ww,pyw;
double rhp[][]= {{0,0,0,0,0,0},{0,0,0,0,0,0},{0,0,0,0,0,0},{0,0,0,0,0,0},
{0,0,0,0,0,0},{0,0,0,0,0,0},{0,0,0,0,0,0},{0,0,0,0,0,0}};
double rpp[][]= {{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}};
double mmp[][]={{0,0,0,0},{0,0,0,0},{0,0,0,0}};
double teqq[][]={{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0}};
double hevec[][]=new double[3][4];
for (int jou=0; jou <=2; jou++) {
hevec[0][jou]=nux[1][jou]-nux[0][jou]; // hevec[0] = O nuc - H0 nuc vector
hevec[1][jou]=nux[2][jou]-nux[0][jou]; // hevec[1] = O nuc - H1 nuc vector
}
// hevec[2][0-2] = vector perpendicular to H-O-H plane = cross product 0 x 1
hevec[2][0]=hevec[0][1]*hevec[1][2]-hevec[0][2]*hevec[1][1];
hevec[2][1]=hevec[0][2]*hevec[1][0]-hevec[0][0]*hevec[1][2];
hevec[2][2]=hevec[0][0]*hevec[1][1]-hevec[0][1]*hevec[1][0];
// hevec[][3]=length of each vector
hevec[2][3]=Math.sqrt(hevec[2][0]*hevec[2][0]+hevec[2][1]*hevec[2][1]+hevec[2][2]*hevec[2][2]);
hevec[0][3]=Math.sqrt(hevec[0][0]*hevec[0][0]+hevec[0][1]*hevec[0][1]+hevec[0][2]*hevec[0][2]);
hevec[1][3]=Math.sqrt(hevec[1][0]*hevec[1][0]+hevec[1][1]*hevec[1][1]+hevec[1][2]*hevec[1][2]);
double noxx[]={0,0,0}; // noxx[0-2]=center of six valence electrons (0-5)
for (int yp=0; yp <=5; yp++) {
for (int jou=0; jou <=2; jou++) {
noxx[jou]=noxx[jou]+hpr[yp][jou];
}}
for (int jou=0; jou <=2; jou++) {
noxx[jou]=noxx[jou]/6.0;
}
// calculate hpr2[6,7][0-2]= symmetric positions of ele 6,7 with respect to H-O-H plane
for (int yp=6; yp <=7; yp++) {
gx=(hevec[2][0]*(hpr[yp][0]-nux[0][0])+hevec[2][1]*(hpr[yp][1]-nux[0][1])
+hevec[2][2]*(hpr[yp][2]-nux[0][2]))/hevec[2][3];
for (int jou=0; jou <=2; jou++) {
hpr2[yp][jou]=hpr[yp][jou]-(gx*2.0*hevec[2][jou])/hevec[2][3];
} }
toav=0.0; // toav=total potential energy (eV)
// kro= distance between two electrons (6 and 7)
kro=Math.sqrt((hpr[6][0]-hpr[7][0])*(hpr[6][0]-hpr[7][0])+
(hpr[6][1]-hpr[7][1])*(hpr[6][1]-hpr[7][1])+
(hpr[6][2]-hpr[7][2])*(hpr[6][2]-hpr[7][2]));
if (kro==0) {kro=5000.0;}
// pot=potential energy between electrons (6 and 7)
pot=(elc*elc*6.241509e18)/(4*pai*epsi*kro*1.0e-14);
// rhp[][3]=potential energy of each electron (eV)
rhp[6][3]=rhp[6][3]+pot/2; rhp[7][3]=rhp[7][3]+pot/2;
toav=toav+pot;
for (int jou=0; jou <=2; jou++) { //ggx=force between two electrons (6 and 7)
ggx=(suh*(hpr[6][jou]-hpr[7][jou]))/(kro*kro*kro);
// rhp[][0-2]=force component of each electron
rhp[6][jou]=rhp[6][jou]+ggx; rhp[7][jou]=rhp[7][jou]-ggx;
}
for (int yp=0; yp <=5; yp++) { // interaction between electrons (0-5) and (6,7)
for (int kj=6; kj <=7; kj++) {
kro=Math.sqrt((hpr[yp][0]-hpr[kj][0])*(hpr[yp][0]-hpr[kj][0])+
(hpr[yp][1]-hpr[kj][1])*(hpr[yp][1]-hpr[kj][1])+
(hpr[yp][2]-hpr[kj][2])*(hpr[yp][2]-hpr[kj][2]));
if (kro==0) {kro=5000.0;}
// kro2=distance between electrons (0-5) and symmetric positions of ele 6,7
kro2=Math.sqrt((hpr[yp][0]-hpr2[kj][0])*(hpr[yp][0]-hpr2[kj][0])+
(hpr[yp][1]-hpr2[kj][1])*(hpr[yp][1]-hpr2[kj][1])+
(hpr[yp][2]-hpr2[kj][2])*(hpr[yp][2]-hpr2[kj][2]));
if (kro2==0) {kro2=5000.0;}
pot=(elc*elc*6.241509e18)/(4*pai*epsi*kro*1.0e-14);
pot2=(elc*elc*6.241509e18)/(4*pai*epsi*kro2*1.0e-14);
rhp[yp][3]=rhp[yp][3]+pot/4.0+pot2/4.0; rhp[kj][3]=rhp[kj][3]+pot/2.0;
toav=toav+pot;
for (int jou=0; jou <=2; jou++) {
ggx=(suh*(hpr[yp][jou]-hpr[kj][jou]))/(kro*kro*kro);
ggy=(suh*(hpr[yp][jou]-hpr2[kj][jou]))/(kro2*kro2*kro2);
rhp[yp][jou]=rhp[yp][jou]+ggx*0.5+ggy*0.5; rhp[kj][jou]=rhp[kj][jou]-ggx;
// rpp[0-5][0-2]= the sum of force acting on ele 0-5 from 3 nuclei and ele 6,7
rpp[yp][jou]=rpp[yp][jou]+ggx*0.5+ggy*0.5;
}
}}
double ppot;
for (int yp=0; yp <=5; yp++) {
for (int kj=0; kj <=5; kj++) {
if (yp < kj ) {
kro=Math.sqrt((hpr[yp][0]-hpr[kj][0])*(hpr[yp][0]-hpr[kj][0])+
(hpr[yp][1]-hpr[kj][1])*(hpr[yp][1]-hpr[kj][1])+
(hpr[yp][2]-hpr[kj][2])*(hpr[yp][2]-hpr[kj][2]));
if (kro==0) {kro=5000.0;}
ppot=(elc*elc*6.241509e18)/(4*pai*epsi*kro*1.0e-14);
toav=toav+ppot;
rhp[yp][3]=rhp[yp][3]+ppot/2.0; rhp[kj][3]=rhp[kj][3]+ppot/2.0;
// teqq[0-5][3] = potential energy only in oxygen atom
teqq[yp][3]=teqq[yp][3]+ppot/2; teqq[kj][3]=teqq[kj][3]+ppot/2;
for (int jou=0; jou <=2; jou++) {
ggx=(suh*(hpr[yp][jou]-hpr[kj][jou]))/(kro*kro*kro);
rhp[yp][jou]=rhp[yp][jou]+ggx; rhp[kj][jou]=rhp[kj][jou]-ggx;
// teqq[0-5][0-2] = force component only in oxygen atom
teqq[yp][jou]=teqq[yp][jou]+ggx; teqq[kj][jou]=teqq[kj][jou]-ggx;
}
}}}
for (int rv=0; rv <=5; rv++) { // interaction between noxx(center of ele 0-5) and ele0-5
kro=Math.sqrt((hpr[rv][0]-noxx[0])*(hpr[rv][0]-noxx[0])+
(hpr[rv][1]-noxx[1])*(hpr[rv][1]-noxx[1])+
(hpr[rv][2]-noxx[2])*(hpr[rv][2]-noxx[2]));
if (kro == 0) {kro=5000.0;}
ppot=(elc*elc*den*6.241509e18)/(4*pai*epsi*kro*1.0e-14);
// teqq[0-5][3] = potential energy only in oxygen atom
teqq[rv][3]=teqq[rv][3]-ppot;
for (int jou=0; jou <=2; jou++) {
ggx=(suh*den*(hpr[rv][jou]-noxx[jou]))/(kro*kro*kro);
// teqq[0-5][0-2] = force component only in oxygen atom
teqq[rv][jou]=teqq[rv][jou]-ggx;
}}
// interaction between electrons and nuclei
for (int rv=0; rv <=7; rv++) {
// -------------------- O nucleus (nux[0][0-2])
kro=Math.sqrt((hpr[rv][0]-nux[0][0])*(hpr[rv][0]-nux[0][0])+
(hpr[rv][1]-nux[0][1])*(hpr[rv][1]-nux[0][1])+
(hpr[rv][2]-nux[0][2])*(hpr[rv][2]-nux[0][2])); kro2=1000.0;
if (kro == 0) {kro=5000.0;}
// kro2 = distance between O nucleus and symmetric positions of ele 6,7
if (rv > 5) {kro2=Math.sqrt((hpr2[rv][0]-nux[0][0])*(hpr2[rv][0]-nux[0][0])+
(hpr2[rv][1]-nux[0][1])*(hpr2[rv][1]-nux[0][1])+
(hpr2[rv][2]-nux[0][2])*(hpr2[rv][2]-nux[0][2]));
if (kro2 == 0) {kro2=5000.0;}
}
pot =-den/kro; // show distance (nuc) between electron 0-5 and O nucleus
if (rv < 6) {ex=(int)(kro); elp[rv][3].setText("nuc "+Integer.toString(ex));}
for (int kj=0; kj <=2; kj++) {
ggx=(suh*den*(hpr[rv][kj]-nux[0][kj]))/(kro*kro*kro); ggy=ggx;
if (rv > 5) {ggy=(suh*den*(hpr2[rv][kj]-nux[0][kj]))/(kro2*kro2*kro2);}
// mmp[0][0-2] = force component acting on O nucleus
mmp[0][kj]=mmp[0][kj]+ggx/2.0+ggy/2.0; rhp[rv][kj]=rhp[rv][kj]-ggx;
// rpp[0-5][0-2]=the sum of force acting on ele 0-5 from 3 nuclei and ele 6,7
if (rv < 6) {rpp[rv][kj]=rpp[rv][kj]-ggx;}
}
// ---------------------------- H0 nucleus (nux[1][0-2])
krr=Math.sqrt((hpr[rv][0]-nux[1][0])*(hpr[rv][0]-nux[1][0])+
(hpr[rv][1]-nux[1][1])*(hpr[rv][1]-nux[1][1])+
(hpr[rv][2]-nux[1][2])*(hpr[rv][2]-nux[1][2]));
if (krr ==0) {krr=5000.0;}
if (rv > 5) {krr2=Math.sqrt((hpr2[rv][0]-nux[1][0])*(hpr2[rv][0]-nux[1][0])+
(hpr2[rv][1]-nux[1][1])*(hpr2[rv][1]-nux[1][1])+
(hpr2[rv][2]-nux[1][2])*(hpr2[rv][2]-nux[1][2]));
if (krr2 == 0) {krr2=5000.0;}
}
pota=-1.0/krr; // show distance (nuc) between electron 6 and H0 nucleus
if (rv==6) {ex=(int)(krr); elp[rv][3].setText("nuc "+Integer.toString(ex));}
for (int kj=0; kj <=2; kj++) {
ggx=(suh*(hpr[rv][kj]-nux[1][kj]))/(krr*krr*krr); ggy=ggx;
if (rv > 5) {ggy=(suh*(hpr2[rv][kj]-nux[1][kj]))/(krr2*krr2*krr2);}
// mmp[1][0-2] = force component acting on H0 nucleus
mmp[1][kj]=mmp[1][kj]+ggx/2.0+ggy/2.0; rhp[rv][kj]=rhp[rv][kj]-ggx;
// rpp[0-5][0-2]=the sum of force acting on ele 0-5 from 3 nuclei and ele 6,7
if (rv < 6) {rpp[rv][kj]=rpp[rv][kj]-ggx;}
}
// ---------------------------- H1 nucleus (nux[2][0-2])
krk=Math.sqrt((hpr[rv][0]-nux[2][0])*(hpr[rv][0]-nux[2][0])+
(hpr[rv][1]-nux[2][1])*(hpr[rv][1]-nux[2][1])+
(hpr[rv][2]-nux[2][2])*(hpr[rv][2]-nux[2][2]));
if (krk ==0) {krk=5000.0;}
if (rv > 5) {krk2=Math.sqrt((hpr2[rv][0]-nux[2][0])*(hpr2[rv][0]-nux[2][0])+
(hpr2[rv][1]-nux[2][1])*(hpr2[rv][1]-nux[2][1])+
(hpr2[rv][2]-nux[2][2])*(hpr2[rv][2]-nux[2][2]));
if (krk2 == 0) {krk2=5000.0;}
}
potb=-1.0/krk; // show distance (nuc) between electron 7 and H1 nucleus
if (rv==7) {ex=(int)(krk); elp[rv][3].setText("nuc "+Integer.toString(ex));}
for (int kj=0; kj <=2; kj++) {
ggx=(suh*(hpr[rv][kj]-nux[2][kj]))/(krk*krk*krk); ggy=ggx;
if (rv > 5) {ggy=(suh*(hpr2[rv][kj]-nux[2][kj]))/(krk2*krk2*krk2);}
// mmp[2][0-2] = force component acting on H1 nucleus
mmp[2][kj]=mmp[2][kj]+ggx/2.0+ggy/2.0; rhp[rv][kj]=rhp[rv][kj]-ggx;
// rpp[0-5][0-2]=the sum of force acting on ele 0-5 from 3 nuclei and ele 6,7
if (rv < 6) {rpp[rv][kj]=rpp[rv][kj]-ggx;}
}
// ttav=potential energy between one electron and three nuclei
ttav= (elc*elc*6.241509e18*(pot+pota+potb))/(4*pai*epsi*1.0e-14);
rhp[rv][3]=rhp[rv][3]+ttav; toav=toav+ttav;
// rhp[][5]= potential energy between electron and other nuclei
if (rv < 6) {rhp[rv][5]=(elc*elc*6.241509e18*(pota+potb))/(4*pai*epsi*1.0e-14); }
if (rv == 6) {rhp[rv][5]=(elc*elc*6.241509e18*(pot+potb))/(4*pai*epsi*1.0e-14);}
if (rv == 7) {rhp[rv][5]=(elc*elc*6.241509e18*(pot+pota))/(4*pai*epsi*1.0e-14);}
}
// kro=distance between O and H0 nuclei
kro=Math.sqrt((nux[0][0]-nux[1][0])*(nux[0][0]-nux[1][0])+
(nux[0][1]-nux[1][1])*(nux[0][1]-nux[1][1])+
(nux[0][2]-nux[1][2])*(nux[0][2]-nux[1][2]));
// forces between O and H0 nuclei
ggx=(suh*den*(nux[1][0]-nux[0][0]))/(kro*kro*kro); ggy=(suh*den*(nux[1][1]-nux[0][1]))/(kro*kro*kro);
ggz=(suh*den*(nux[1][2]-nux[0][2]))/(kro*kro*kro);
mmp[0][0]=mmp[0][0]-ggx; mmp[0][1]=mmp[0][1]-ggy; mmp[0][2]=mmp[0][2]-ggz;
mmp[1][0]=mmp[1][0]+ggx; mmp[1][1]=mmp[1][1]+ggy; mmp[1][2]=mmp[1][2]+ggz;
// krr=distance between O and H1 nuclei
krr=Math.sqrt((nux[0][0]-nux[2][0])*(nux[0][0]-nux[2][0])+
(nux[0][1]-nux[2][1])*(nux[0][1]-nux[2][1])+
(nux[0][2]-nux[2][2])*(nux[0][2]-nux[2][2]));
// forces between O and H1 nuclei
ggx=(suh*den*(nux[2][0]-nux[0][0]))/(krr*krr*krr); ggy=(suh*den*(nux[2][1]-nux[0][1]))/(krr*krr*krr);
ggz=(suh*den*(nux[2][2]-nux[0][2]))/(krr*krr*krr);
mmp[0][0]=mmp[0][0]-ggx; mmp[0][1]=mmp[0][1]-ggy; mmp[0][2]=mmp[0][2]-ggz;
mmp[2][0]=mmp[2][0]+ggx; mmp[2][1]=mmp[2][1]+ggy; mmp[2][2]=mmp[2][2]+ggz;
krk=nux[2][1]-nux[1][1];
// forces between H0 and H1 nuclei
ggy=(suh)/(krk*krk);
mmp[1][1]=mmp[1][1]-ggy; mmp[2][1]=mmp[2][1]+ggy;
// pot=repulsive potential energy among three nuclei (eV)
pot=(elc*elc*6.241509e18*((den/kro)+(den/krr)+(1.0/krk)))/(4*pai*epsi*1.0e-14);
toav=toav+pot;
ex=(int)(100*toav); ggx=ex/100.0;
impho.setText("tV "+Double.toString(ggx)); // show total V to two decimal places
// distribute repulsive V among nuclei to each electron based on rhp[][5]
double hiwa=0.0;
for (int rv=0; rv <=7; rv++) { hiwa=hiwa+rhp[rv][5]; }
for (int rv=0; rv <=7; rv++) {
rhp[rv][3]=rhp[rv][3]+(pot*rhp[rv][5])/hiwa;
ex=(int)(100*rhp[rv][3]); ggx=ex/100.0;
elp[rv][4].setText("V "+Double.toString(ggx)); // show each electron's V (=rhp[][3])
}
gx=0.0; // gx=sum of each potential energy
for (int rv=0; rv <=7; rv++) { gx=gx+rhp[rv][3]; }
gy=-toav*0.5; // gy=total kinetic energy
// distribute kinetic energy to each electron based on rhp[][3]
for (int rv=0; rv <=7; rv++) {
gz=(gy*rhp[rv][3])/gx; rhp[rv][4]=gz; // rhp[][4]=each electron's kinetic energy
ex=(int)(100*gz); gz=ex/100.0; elp[rv][5].setText("T "+Double.toString(gz));
}
for (int rv=0; rv <=5; rv++) { // show force component acting on electron 0-5
gx=Math.sqrt(rpp[rv][0]*rpp[rv][0]+rpp[rv][1]*rpp[rv][1]+rpp[rv][2]*rpp[rv][2]);
// gy=force component (CF) in the direction of rpp = inner product of rhp and rpp
gy=(rhp[rv][0]*rpp[rv][0]+rhp[rv][1]*rpp[rv][1]+rhp[rv][2]*rpp[rv][2])/gx;
ex=(int)(1000*gy); ww="CF ";
elp[rv][6].setText(ww+Integer.toString(ex));
for (int jou=0; jou <=2; jou++) {
gz=rhp[rv][jou]-(gy*rpp[rv][jou])/gx;
ex=(int)(1000*gz); // show force component other than CF
elp[rv][jou+7].setText(Integer.toString(ex));
}
}
for (int rv=6; rv <=7; rv++) { // show force component acting on electron 6 and 7
// gy= force component(CF) of ele6,7 peperndicular to H-O-H plane = inner product of rhp and hevec[2]
gy=(rhp[rv][0]*hevec[2][0]+rhp[rv][1]*hevec[2][1]+rhp[rv][2]*hevec[2][2])/hevec[2][3];
ex=(int)(1000*gy); ww="CF ";
elp[rv][6].setText(ww+Integer.toString(ex));
for (int jou=0; jou <=2; jou++) {
gz=rhp[rv][jou]-(gy*hevec[2][jou])/hevec[2][3];
ex=(int)(1000*gz); // show force component other than CF
elp[rv][jou+7].setText(Integer.toString(ex));
}
}
// mmp[0][3] = force component of O nucleus toward H0 nucleus
mmp[0][3]=(mmp[0][0]*hevec[0][0]+mmp[0][1]*hevec[0][1]+mmp[0][2]*hevec[0][2])/hevec[0][3];
// mmp[1][3] = force component of H0 nucleus toward O nucleus
mmp[1][3]=-(mmp[1][0]*hevec[0][0]+mmp[1][1]*hevec[0][1]+mmp[1][2]*hevec[0][2])/hevec[0][3];
// mmp[2][3] = force component of H1 nucleus toward O nucleus
mmp[2][3]=-(mmp[2][0]*hevec[1][0]+mmp[2][1]*hevec[1][1]+mmp[2][2]*hevec[1][2])/hevec[1][3];
for (int rv=0; rv <=2; rv++) {
for (int jou=0; jou <=3; jou++) { // show mmp[0-2][0-3] = force acting on each nucleus
ex=(int)(1000*mmp[rv][jou]); ww=" ";
if (jou==0) {ww="FX=";}
if (jou==1) {ww="FY=";}
if (jou==2) {ww="FZ=";}
if (jou==3 && rv==0) {ww="H0 ";} if (jou==3 && rv > 0) {ww="0n ";}
mmpho[rv][jou].setText(ww+Integer.toString(ex));
}}
for (int rv=0; rv <=7; rv++) { // show de Broglie wave of each electron
gz=Math.sqrt(rhp[rv][0]*rhp[rv][0]+rhp[rv][1]*rhp[rv][1]+rhp[rv][2]*rhp[rv][2]);
// electrons 0-5 use forces (=teqq[0-5][0-2]) only in oxygen atom
if (rv < 6) {gz=Math.sqrt(teqq[rv][0]*teqq[rv][0]+teqq[rv][1]*teqq[rv][1]+teqq[rv][2]*teqq[rv][2]);}
gy=(gz*elc*elc)/(4*pai*epsi*suh*1.0e-28); // gy=force (N)
gx=Math.sqrt((2*rhp[rv][4]*1.602177e-19)/me); // gx=velocity (m/s) from kinetic energy
// electrons 0-5 use potential V (and T) only in oxygen atom
if (rv < 6) {gx=Math.sqrt((-teqq[rv][3]*1.602177e-19)/me);}
ggx=(me*gx*gx)/gy; // ggx= "tenporary" radius (m)
ggy=(2*pai*ggx*me*gx)/h; // ggy (wn) = number of de Broglie's waves contained in one orbit
ex=(int)(ggy*1000); ggy=ex/1000.0; // show wn to three decimal places
elp[rv][10].setText("wn "+Double.toString(ggy)); }
// --------------------- show picture
int nmx[][]=new int[3][3]; int hpk[][]=new int[8][4];
for (int yp=0; yp <=2; yp++) {
for (int kj=0; kj <=2; kj++) { // change MM to pixel in nuclei
nmx[yp][kj]=(int)(nux[yp][kj]/71.428);
}}
for (int yp=0; yp <=7; yp++) {
for (int kj=0; kj <=2; kj++) { // change MM to pixel in electrons
hpk[yp][kj]=(int)(hpr[yp][kj]/71.428);
}}
g.clearRect(9,299,1170,699);
g.setColor(Color.cyan); g.drawLine(375,310,375,660); g.drawLine(735,310,735,660);
g.setColor(Color.lightGray); // show three nuclei
g.fillOval(nmx[0][0]+10,650-nmx[0][1],20,20);g.fillOval(370+nmx[0][0],650-nmx[0][2],20,20);
g.fillOval(730+nmx[0][1],650-nmx[0][2],20,20);
g.fillOval(13+nmx[1][0],653-nmx[1][1],14,14);g.fillOval(373+nmx[1][0],653-nmx[1][2],14,14);
g.fillOval(733+nmx[1][1],653-nmx[1][2],14,14);
g.fillOval(13+nmx[2][0],653-nmx[2][1],14,14);g.fillOval(373+nmx[2][0],653-nmx[2][2],14,14);
g.fillOval(733+nmx[2][1],653-nmx[2][2],14,14);
g.setColor(Color.white); // show electron 0 (particle)
g.fillOval(hpk[0][0]+13,653-hpk[0][1],14,14);
g.fillOval(hpk[0][0]+373,653-hpk[0][2],14,14);
g.fillOval(hpk[0][1]+733,653-hpk[0][2],14,14);
// show electron 1
g.fillOval(hpk[1][0]+13,653-hpk[1][1],14,14);
g.fillOval(hpk[1][0]+373,653-hpk[1][2],14,14);
g.fillOval(hpk[1][1]+733,653-hpk[1][2],14,14);
g.setColor(Color.red); // show electron 2
g.fillOval(hpk[2][0]+13,653-hpk[2][1],14,14);
g.fillOval(hpk[2][0]+373,653-hpk[2][2],14,14);
g.fillOval(hpk[2][1]+733,653-hpk[2][2],14,14);
// show electron 3
g.fillOval(hpk[3][0]+13,653-hpk[3][1],14,14);
g.fillOval(hpk[3][0]+373,653-hpk[3][2],14,14);
g.fillOval(hpk[3][1]+733,653-hpk[3][2],14,14);
// show electron 4
g.setColor(Color.green);
g.fillOval(hpk[4][0]+13,653-hpk[4][1],14,14);
g.fillOval(hpk[4][0]+373,653-hpk[4][2],14,14);
g.fillOval(hpk[4][1]+733,653-hpk[4][2],14,14);
// show electron 5
g.fillOval(hpk[5][0]+13,653-hpk[5][1],14,14);
g.fillOval(hpk[5][0]+373,653-hpk[5][2],14,14);
g.fillOval(hpk[5][1]+733,653-hpk[5][2],14,14);
// show electron 6
g.setColor(Color.pink);
g.fillOval(hpk[6][0]+13,653-hpk[6][1],14,14);
g.fillOval(hpk[6][0]+373,653-hpk[6][2],14,14);
g.fillOval(hpk[6][1]+733,653-hpk[6][2],14,14);
// show electron 7
g.fillOval(hpk[7][0]+13,653-hpk[7][1],14,14);
g.fillOval(hpk[7][0]+373,653-hpk[7][2],14,14);
g.fillOval(hpk[7][1]+733,653-hpk[7][2],14,14);
for (int rw=0; rw <=7; rw++) { // show each electron's number
g.setColor(Color.blue);
g.drawString(Integer.toString(rw),hpk[rw][0]+17,665-hpk[rw][1]);
g.drawString(Integer.toString(rw),hpk[rw][0]+377,665-hpk[rw][2] );
g.drawString(Integer.toString(rw),hpk[rw][1]+737,665-hpk[rw][2] );
}
}
}