SUBROUTINE calc_pot2 (nrho, ntheta, nchi, rho_val, theta_val, chi_val, vpot, morph, morph_val) USE Numeric_Kinds_Module USE Numbers_Module USE FileUnits_Module USE PES_MODULE USE Masses_Module USE InputFile_Module IMPLICIT NONE LOGICAL :: const_rho, const_chi, const_theta, vmin_rho, vmin_theta, There LOGICAL :: vmin_chi, morph LOGICAL, Parameter :: Prt_Masses=.False. INTEGER :: nrho, nchi, ntheta, irho, itheta, ichi, icase, nkeyfr, kfrcnt, I, IERR REAL(Kind=WP_Kind) :: rho_min, rho_max, rho, chi,theta, upot, drho REAL(Kind=WP_Kind) :: dchi, theta_min, theta_max, chi_min, chi_max REAL(Kind=WP_Kind) :: vmin, radtodeg, morph_val REAL(Kind=WP_Kind) :: rho_val (nrho), theta_val (ntheta), chi_val (nchi) REAL(Kind=WP_Kind) :: vpot (nrho, ntheta, nchi) REAL(Kind=WP_Kind), ALLOCATABLE:: temp(:,:) NAMELIST/pot2/ rho_min, rho_max, const_rho, rho, vmin_rho, const_chi, chi, vmin_chi, const_theta, theta, vmin_theta ALLOCATE(temp(ntheta,ntheta)) WRITE(Out_Unit, * ) WRITE(Out_Unit, * ) 'Begin(Calc_Pot2) Calculate Potential' ! Atomic Masses INQUIRE (File=InputDIR(1:LEN(TRIM(InputDIR)))//TRIM(InputFile), exist = there) IF(there)THEN REWIND In_Unit OPEN(Unit=In_Unit,File=TRIM(InputDIR)//TRIM(InputFile)) READ(In_Unit, NML=AtomicMasses, IOSTAT=IERR) ! Atomic Masses IF(IERR==-1)THEN IERR=0 REWIND In_Unit READ(In_Unit, NML=Atoms, IOSTAT=IERR) IF(Prt_Masses)THEN WRITE(Out_Unit,NML=Atoms) WRITE(Msg_Unit,NML=Atoms) WRITE(Out_Unit,*) ENDIF DO I=1,3 CALL AtomicWeights(AtomicSymbol(I), MassNumber(I), AtomicNumber(I), Mass(I), Abundance(I), AtomicWeight(I), Notes(I)) IF(Prt_Masses)THEN WRITE(Out_Unit, '(A,I3)')'Atomic Number = ', AtomicNumber(I) WRITE(Out_Unit, '(A,A3)')'Atomic Symbol = ', TRIM(AtomicSymbol(I)) WRITE(Out_Unit, '(A,I3)')'Mass Number = ' , MassNumber(I) WRITE(Out_Unit, '(A,1PES23.15)')'Relative Atomic Mass = ' , Mass(I) WRITE(Out_Unit, '(A,1PES23.15)')'Isotopic Composition = ' , Abundance(I) WRITE(Out_Unit, '(A,1PES23.15)')'Standard Atomic Weight = ', AtomicWeight(I) WRITE(Out_Unit, '(A,A)')'Notes = ', TRIM(Notes(I)) WRITE(Out_Unit,*) ENDIF ENDDO ELSEIF (IERR/=0)THEN WRITE(Msg_Unit,*)' Missing a necessary namelist: Atoms or AtomicMasses' ENDIF WRITE(Out_Unit,NML=AtomicMasses) CLOSE(In_Unit) ELSE WRITE(Out_Unit, * ) 'Data file must exist: ', InputDIR(1:LEN(TRIM(InputDIR)))//InputFile WRITE(Msg_Unit, * ) 'Data file must exist: ', InputDIR(1:LEN(TRIM(InputDIR)))//InputFile STOP "Calc_Pot2" ENDIF OPEN(Unit=In_Unit,File=InputDIR(1:LEN(TRIM(InputDIR)))//InputFile) READ(In_Unit,NML=pot2, IOSTAT=IERR) IF(IERR==-1)THEN CALL NameList_Default(Out_Unit, 'Pot2') rho_min=1.d0 rho_max=10.d0 const_rho=.false. rho=9.5d0 vmin_rho=.false. const_theta=.false. theta=Pi/two vmin_theta=.false. const_chi=.false. chi=0.0 vmin_chi=.false. ELSEIF(IERR/=0)THEN WRITE(Out_Unit,*)'ERROR with Namelist Pot2' WRITE(Msg_Unit,*)'ERROR with Namelist Pot2' STOP 'Calc_Pot2: Pot2' ENDIF CLOSE(In_Unit) WRITE(Out_Unit,NML=pot2) CALL constant CALL massfac (.true.) IF(const_rho)THEN nrho = 1 IF(morph)THEN rho = morph_val ENDIF rho_val (1) = rho ELSE drho = (rho_max - rho_min) / (nrho - 1) DO irho = 1, nrho rho_val (irho) = (irho - 1) * drho + rho_min ENDDO ENDIF theta_min = 0.d0 theta_max = 2.d0 * atan (1.d0) IF(const_theta)THEN ntheta = 1 IF(morph)THEN theta = morph_val ENDIF theta_val (1) = theta ELSE !-------------------------------------- !set new theta values !-------------------------------------- CALL DVR_kinetic(ntheta,temp,theta_val,One) ENDIF chi_min = 0.d0 chi_max = 8.d0 * atan (1.d0) IF(const_chi)THEN nchi = 1 IF(morph)THEN chi = morph_val ENDIF chi_val (1) = chi ELSE dchi = (chi_max - chi_min) / (nchi - 1) DO ichi = 1, nchi chi_val (ichi) = (ichi - 1) * dchi ENDDO ENDIF kfrcnt = 0 DO ichi = 1, nchi chi = chi_val (ichi) DO itheta = 1, ntheta theta = theta_val (itheta) DO irho = 1, nrho rho = rho_val (irho) nkeyfr = nchi * ntheta * nrho kfrcnt = kfrcnt + 1 CALL potaph2 (rho, theta, chi, upot, nkeyfr, kfrcnt) vpot (irho, itheta, ichi) = upot ENDDO ENDDO ENDDO IF(kfrcnt/=nkeyfr)THEN STOP 'kfrcnt_calc_pot.f' ENDIF IF(vmin_rho)THEN DO itheta = 1, ntheta DO ichi = 1, nchi vmin = Ten**(RANGE(One)-10) DO irho = 1, nrho IF(vpot (irho, itheta, ichi)