!! MODULDE provides spin-independent separable rank one interaction in S-wave 

MODULE yamaguchi
 USE precision
 USE constants
 IMPLICIT NONE
 PRIVATE
 REAL(dpreal) :: lambda1s0,beta1s0,lambda3s1,beta3s1
 REAL(dpreal) :: mfaktor=hbarc/41.47
 LOGICAL :: first=.true.
 
 PUBLIC potyglsj
 
 CONTAINS
 
 SUBROUTINE readpara
  IMPLICIT NONE
  ! read parameter sets 
  OPEN(UNIT=75,FILE='yamaguchi-para.dat',FORM='FORMATTED',STATUS='OLD')
  READ(75,*,END=100,ERR=100)
  READ(75,*,END=100,ERR=100) lambda1s0,beta1s0,lambda3s1,beta3s1                
  CLOSE(75)
  
  WRITE(*,*) 'Yamaguchi separable interaction'
  WRITE(*,*) 'lambda(1S0),beta(1S0) = ',lambda1s0,beta1s0
  WRITE(*,*) 'lambda(3S1),beta(3S1) = ',lambda3s1,beta3s1 
  
  RETURN
  
  100 CONTINUE
  WRITE(*,*) 'problem reading yamaguchi-para.dat'
  STOP
  
 END SUBROUTINE readpara
 
 SUBROUTINE potyglsj(P,PS,L,IS,JJ,V00,V02,V20,V22)
  IMPLICIT NONE
  REAL(dpreal) :: V00,V02,V20,V22,P,PS
  INTEGER L,IS,JJ
  IF(first) THEN
   CALL readpara
   first=.false.
  END IF
  
  IF(JJ.EQ.0 .AND. IS.EQ.0 .AND. L.EQ.0)  THEN ! 1s0
    V00=lambda1s0/((P**2+beta1s0**2)*(PS**2+beta1s0**2))
  ELSE IF(JJ.EQ.1 .AND. IS.EQ.1 .AND. L.EQ.0) THEN ! 3s1 
    V00=lambda3s1/((P**2+beta3s1**2)*(PS**2+beta3s1**2))
  ELSE
   V00=0.0_dpreal 
  END IF
  
  V00=mfaktor*P*PS*V00
  V02=0.0_dpreal
  V20=0.0_dpreal
  V22=0.0_dpreal
 END SUBROUTINE potyglsj 
  
END MODULE yamaguchi