SUBROUTINE pot(r, hstep, v0, v1, v2, jvalue, mu)
USE Numeric_Kinds_Module
USE FileUnits_OneDim_Module
USE FileUnits_Common_Module
USE OneDim_Parms_Module
USE DiatomicPot_MODULE
USE Numbers_Module
IMPLICIT NONE
INTEGER jvalue
REAL(KIND=WP_Kind) r1(3), V0_ph, V0_mh
REAL(KIND=WP_Kind), PARAMETER :: HydrogenMass= 1.00782503207_WP_Kind
REAL(KIND=WP_Kind), PARAMETER :: DeuteriumMass=2.0141017778_WP_Kind
REAL(KIND=WP_Kind), PARAMETER :: TritiumMass=  3.0160492777_WP_Kind
REAL(KIND=WP_Kind), PARAMETER :: H4Mass=       4.02781_WP_Kind
REAL(KIND=WP_Kind), PARAMETER :: H5Mass=       5.03531_WP_Kind
REAL(KIND=WP_Kind), PARAMETER :: H6Mass=       6.04494_WP_Kind
REAL(KIND=WP_Kind), PARAMETER :: H7Mass=       7.05275_WP_Kind
REAL(KIND=WP_Kind) r, hstep, h, v0, v1, v2, mu, cent0, cent1, cent2, VNaCs

IF(HydrogenMass/=MassA.or.HydrogenMass/=MassB)THEN
   WRITE(Msg_Unit,*)"HydrogenMass/=MassA.or.HydrogenMass/=MassB"
   WRITE(Msg_Unit,*)HydrogenMass, MassA, HydrogenMass, MassB
   WRITE(Out_Unit,*)"HydrogenMass/=MassA.or.HydrogenMass/=MassB"
   WRITE(Out_Unit,*)HydrogenMass, MassA, HydrogenMass, MassB
   STOP "Pot - Masses"
ENDIF

IF(DiatomicPot/="H2_PKH3")THEN
   WRITE(Msg_Unit,*)"DiatomicPot/=H2_PKH3"
   WRITE(Msg_Unit,*)DiatomicPot, "H2_PKH3"
   WRITE(Out_Unit,*)"DiatomicPot/=H2_PKH3"
   WRITE(Out_Unit,*)DiatomicPot, "H2_PKH3"
   STOP "Pot - DiatomicPot"
ENDIF

mu=1822.888506d0*HydrogenMass*HydrogenMass/(HydrogenMass+HydrogenMass)
r1(1)=r
r1(2)=1.d+30
r1(3)=1.d+30
h=1.0d-3
CALL H3PK(V0,r1)
r1(1)=r+h
CALL H3PK(V0_ph,r1)
r1(1)=r-h
CALL H3PK(V0_mh,r1)
v1=(V0_ph-V0_mh)/h
v2=(V0_ph-Two*V0+V0_mh)/(Two*h**2)

CALL CenTerm(r, jvalue, mu, cent0, cent1, cent2)
v0=v0+cent0
v1=v1+cent1
v2=v2+cent2

RETURN
END SUBROUTINE pot