ABA or ABC type molecule-project introduction, not checked values

'Needed array buildup for A B C dynamics

rAB[2]

rBC[2]

xA[2]

mA

mB

xB[2]

xC[2]

8.154417198599252`*^-8

'one configuration may be as follows

Do[{
xA[n+1]=xA[n]+h * vA[n]+(h^2)* fA[n]/(2 mA),
xB[n+1]=xB[n]+h * vB[n]+(h^2)* fB[n]/(2 mB),
xC[n+1]=xC[n]+h * vC[n]+(h^2)* fC[n]/(2 mC),
rAB[n+1]=xB[n+1]-xA[n+1],
rBC[n+1]=xC[n+1]-xB[n+1],
fA[n+1]=De (-2 (E^(-2 (rAB[n+1]-re) ββ)) ββ+2 (E^(-(rAB[n+1]-re) ββ)) ββ) ,
fB[n+1]=-(De (-2 (E^(-2 (rAB[n+1]-re) ββ)) ββ+2 (E^(-(rAB[n+1]-re) ββ)) ββ))+(De (-2 (E^(-2 (rBC[n+1]-re) ββ)) ββ+2 (E^(-(rBC[n+1]-re) ββ)) ββ)),
fC[n+1]=-(De (-2 (E^(-2 (rBC[n+1]-re) ββ)) ββ+2 (E^(-(rBC[n+1]-re) ββ)) ββ)),
vA[n+1]=vA[n]+h * (fA[n+1]+fA[n])/(2   mA),
vB[n+1]=vB[n]+h * (fB[n+1]+fB[n])/(2   mB),
vC[n+1]=vC[n]+h * (fC[n+1]+fC[n])/(2   mC)},
{n,1,150,1}]//N

'time interval

h

'possible range of initial conditions\

xB[1]=-0.001 10^(-10)

xA[1]=-1.2751 10^(-10)

xC[1]=+1.2751 10^(-10)

re=1.275  (10^(-10))

mA=1*1836*9.1 10^(-31)

mB=35.5*1836*9.1 10^(-31)

mC=mA

'tentative results and visualizations

vbt=Table[vB[i],{i,1,150}]

vat=Table[vA[i],{i,1,150}]

vct=Table[vC[i],{i,1,150}]

vA[1]

vB[1]

g3=ListPlot[{vct},PlotJoined->True]

g1=ListPlot[{vat},PlotJoined->True]

g2=ListPlot[{vbt},PlotJoined->True]

Show[g1,g2,g3]

ContourPlot[-(De (-2 (E^(-2 (rABx-re) ββ)) ββ+2 (E^(-(rABx-re) ββ)) ββ))+(De (-2 (E^(-2 (rBCx-re) ββ)) ββ+2 (E^(-(rBCx-re) ββ)) ββ)),{rBCx,1. 10^-10,1.76 10^-10},{rABx,1.1 10^-10,1.76 10^-10}]

Plot3D[-(De (-2 (E^(-2 (rABx-re) ββ)) ββ+2 (E^(-(rABx-re) ββ)) ββ))+(De (-2 (E^(-2 (rBCx-re) ββ)) ββ+2 (E^(-(rBCx-re) ββ)) ββ)),{rBCx,1. 10^-10,1.76 10^-10},{rABx,1. 10^-10,1.76 10^-10},PlotStyle->Directive[Yellow,Specularity[White,20],Opacity[0.8]],ExclusionsStyle->{None,Red}]

rABt=Table[rAB[i],{i,1,150}]

rBCt=Table[rBC[i],{i,1,150}]

DensityPlot[{De  (- 2  Exp[-ββ  (rx-re)]+ Exp[-ββ   2 (rx-re)])+De  (- 2  Exp[-ββ  (ry-re)]+ Exp[-ββ   2 (ry-re)])},{rx,1 (10^(-10)),1.77 (10^(-10))},{ry,1.1 (10^(-10)),1.77  (10^(-10))},ColorFunction->Function[{rx,ry,z},Hue[z]]]

DensityPlot[{De  (- 2  Exp[-ββ  (rx-re)]+ Exp[-ββ   2 (rx-re)])+De  (- 2  Exp[-ββ  (ry-re)]+ Exp[-ββ   2 (ry-re)])},{rx,1 (10^(-10)),1.77 (10^(-10))},{ry,1.1 (10^(-10)),1.77  (10^(-10))},ColorFunction->Function[{rx,ry,z},Hue[z]]]

DensityPlot[{De  (- 2  Exp[-ββ  (rx-re)]+ Exp[-ββ   2 (rx-re)])+De  (- 2  Exp[-ββ  (ry-re)]+ Exp[-ββ   2 (ry-re)])},{rx,1 (10^(-10)),1.77 (10^(-10))},{ry,1.1 (10^(-10)),1.77  (10^(-10))}]

ryd=re

Plot3D[{De  (- 2  Exp[-ββ  (rx-re)]+ Exp[-ββ   2 (rx-re)])+De  (- 2  Exp[-ββ  (ry-re)]+ Exp[-ββ   2 (ry-re)])},{rx,1.05 (10^(-10)),2.1 (10^(-10))},{ry,1.07 (10^(-10)),2.28 (10^(-10))}]

ColorFunction->Hue

ContourPlot[-(De (-2 (E^(-2 (rABx-re) ββ)) ββ+2 (E^(-(rABx-re) ββ)) ββ))+(De (-2 (E^(-2 (rBCx-re) ββ)) ββ+2 (E^(-(rBCx-re) ββ)) ββ)),{rBCx,1. 10^-10,1.76 10^-10},{rABx,1.1 10^-10,1.76 10^-10},ColorFunction->Hue]

ColorFunction->Function[{x,y,z},Hue[z]]