!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!System Information!!!!!!!!!!
!!!!LiH2 Reactive Scattering!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 $Title_Labels
   Title1='Li+H2 reactive scattering using the LiH2_PCCP_2015_11732 potential'
   Title2='channel(1)=Li+HH, channel(2)=H+LiH, channel(3)=H+HLi'
   Title3='LiH2_PCCP_2015_11732 PES of Yuan, He, and Chen'
 $END Title_Labels
 
 $PotentialNames
    PES_Name='LiH2_PCCP_2015_11732    ' 
 $end
 
 $FlowOptions
   CallPlotVee    = .False.
   CallBound3D    = .False.
   CallFEM        = .False.
   CallDVR        = .False.
   CallPDAF       = .False.
   CallDelves     = .False.
   CallDiatomic   = .True.
   CalcSFunc      = .True.
   CallABM        = .True.
   CallPropagate  = .True.
   CallAsymptotic = .True.
 $END

 $Atoms  
   AtomicSymbol= 'Li','H','H'
   MassNumber=    7,  1,  1
 $end

 $PlotOpt
   PlotVee=.False.  
 $end
 
 $SectorDef
   DVRBasis        =.False.
   FEMBasis        =.False.
   RhoStart        = 3.23D0
   RhoStartABM_APH = 3.23D0
   RhoStartDelves  = 20.0D0
   RStartJacobi    = 20.0D0
   DeltaRho1       = 0.00D0
   Kase_CalcParm   = 2
   NSectors= 40 !460
 $end

 $regins 
   start        = 3.23d0
   SwitchToVIVS =20.0d0
   finish       =20.0d0
   endaph       =20.0d0
   enddelve     =20.0d0
 $end

 $quantum
  minvib=  0   0  0
  maxvib=  4   4  4
  jmin(0,1,0)=12*0,jmax(0,1)=2   !18  ! Use (maxvib+1)*(jmax or jmin)
  jmin(0,2,0)=12*0,jmax(0,2)=2   !40  ! Use (maxvib+1)*(jmax or jmin)
  jmin(0,3,0)=12*0,jmax(0,3)=2   !40  ! Use (maxvib+1)*(jmax or jmin)
 $end

 $system
  NSfunc=36   !500
 $end

 $Basis_n
  Lam_n=36   !500
  Lamfst=0
  Lamlst=0
 $END
 
 $gauss
   !              HH            LiH            LiH
   npow=           0,             0,             0,
   nglegn=        70,            70,            70,
   nhermt=        70,            70,            70,
   !nlegndre=      18,            40,            40,  !These numbers must be equal to jmax(0,tau)
   intwt=          1,             1,             1,
   !rx=    1.085000D0,    1.075000D0,    1.075000D0,
   !re=    1.401043D0,    3.015247D0,    3.015247D0,
   rx=    1.0850000D0,    1.075000D0,    1.075000D0,
   re=    1.400641D0,    3.015247D0,    3.015247D0,
   zeta=  1.000000D0,    1.000000D0,    1.000000D0,
   delta= 0.010000D0,    0.010000D0,    0.010000D0,
   !weau=  2.0075707D-02, 6.40526285D-3, 6.40526285D-03, 
   !wexeau=5.6068187D-04, 1.055441D-04,  1.055441D-04
   weau=  2.0083744D-02, 6.40526285D-3, 6.40526285D-03,
   wexeau=5.6581021D-04, 1.055441D-04,  1.055441D-04,
   ralpha=3.166000D0,    3.166000D0,    3.166000D0,
   balpha=1.000000D0,    1.000000D0,    1.000000D0,
   calpha=0.920000D0,    0.910000D0,    0.910000D0,
   dalpha=0.030000D0,    0.000000D0,    0.000000D0,
   anharm=0.800000D0,    0.900000D0,    0.900000D0,   
   !old_way=  T
 $end

 $TotEnergy  !low and high energy scattering results
   eV_Input     = .True.
   Efirst_eV    = 0.009D0  ! The total energy is the collision energy for Li+H2 (v=0,j=0)
   DeltaEng_eV  = 0.005D0  ! Emax=Efirst+DeltaEng*(Nenergy-1)=0.001+(.005)*600=3.001                        
   nenergy= 21  !601
 $end

 $momentum  
   jtot=0, 
   parity=0, 
   nsymc=F  F  F, 
   jeven=T, 
   jref=0, 
   megacoup=0, 
   megamax=0
 $end

 $sfplot_ABM
    NTheta_Plot=40
    a=90.d0
    b=0.2d0 
    c=1.d-1
    phi_s=0.5d0
    ithmode=1
    primitive=.false.
    Plot_Sfunc=.False.
 $end

$convrg
   eigmin    = 1.D-5
   ovrerr    = 1.d-6
   ngood     = 15
   nfreq(0)  = 15
   nave1     = 36
   nave2     = 36
 $end

 $convrgxx
   eigmin    = 1.D-2
   ovrerr    = 1.d-3
   ngood     = 500
   nfreq(0)  = 500
   nave1     = 70
   nave2     = 500
 $end

 $intchanl 
   integrat(1)=.true.,
   integrat(2)=.true.,
   integrat(3)=.true.
 $end

 $fbr
   nsfunc=40
   lmax=60
   mmax=120
   nsym=1
 $end

-----------------------------------------------------------------
 Read in the spectroscopic constants to define the experimental vibrational
 rotational states of the isolated diatoms. These constants are obtained
 from Herzberg's book Spectra of Diatomic Molecules and have units of
 inverse centimeters. These are used to insure that one
 has calculated accurate asymptotic wavefunctions.  The calculated energies
 are used and not the 'experimental' energies.
-----------------------------------------------------------------

$spectro
   wecm(1)=   4.416736061650146E+03,    !H2
   wexecm(1)= 1.283607645625079E+02,
   weyecm(1)= 2.359888976108824E+00,
   becm(1)=   6.053858446565800E+01,
   alfecm(1)= 2.851588038673988E+00,
   decm(1)=   4.176459140051390E-02,

   wecm(2)=   1.405776052699711E+03,    !LiH
   wexecm(2)= 2.322485842320376E+01,
   weyecm(2)= 1.895184810192147E-01,
   becm(2)=   7.428559348751051E+00,
   alfecm(2)= 2.151903502000777E-01,
   decm(2)=   6.730958711938912E-04,

   wecm(3)=   1.405776052699711E+03,    !LiH
   wexecm(3)= 2.322485842320376E+01,
   weyecm(3)= 1.895184810192147E-01,
   becm(3)=   7.428559348751051E+00,
   alfecm(3)= 2.151903502000777E-01,
   decm(3)=   6.730958711938912E-04,

 $end


 $spectroxx_old
   wecm(1)=   4406.109551,    !H2
   wexecm(1)=  123.055644,
   weyecm(1)=    1.416041,
   becm(1)=      60.87792,
   alfecm(1)=     3.04684,
   decm(1)=       0.01950,

   wecm(2)=   1405.792698,    !LiH
   wexecm(2)=   23.164268,
   weyecm(2)=    0.174858,
   becm(2)=      7.519308,
   alfecm(2)=    0.216502,
   decm(2)=      0.000000,

   wecm(3)=   1405.792698,    !LiH
   wexecm(3)=   23.164268,
   weyecm(3)=    0.174858,
   becm(3)=      7.519308,
   alfecm(3)=    0.216502,
   decm(3)=      0.000000,
 $end