Acoefficients and new write inelastic probabilities

SRC/AUXILIARy-CODES/distriREACwvp.f

@@ -51,7 +51,7 @@ program distriREACwvp ! version v6 in progress
  real*8,allocatable,dimension(:) :: vibprod,rotprod,vibrot
  integer,allocatable,dimension(:,:) :: noBCstates
  complex*16,allocatable,dimension(:,:,:,:) :: zS2prod
- integer :: noloopreal,ifail
+ integer :: noloopreal,ifail,ncanmax,ican
 
 **>> constants
 
@@ -158,7 +158,8 @@ program distriREACwvp ! version v6 in progress
  j1=jmax
  nv0=nvini
  nv1=nvmax
- allocate(Cvjprod(nv0:nv1,j00:j1,iom0:iom1))
+ ncanmax=(iom1-iom0+1)*nelecmax
+ allocate(Cvjprod(nv0:nv1,j00:j1,ncanmax))
 
  xmasa=(xm2*(xm0+xm1))/xmtot
  xmasaprod=xmasa
@@ -217,15 +218,17 @@ program distriREACwvp ! version v6 in progress
  do ik=1,ntimes
  kcheby=kcheby+1
  read(16,err=2)it,Cvjprod
- do Iom=iom0,iom1
+ ican=0
  do ielec=1,nelecmax
+ do Iom=iom0,iom1
+ ican=ican+1
  do j=j00,j1 
  do iv=nv0,noBCstates(j,ielec)
  if(iv.eq.nvref.and.j.eq.jref.and.iom.eq.iomref
  & .and.kcheby.lt.kminelastic)then
  prodCvj(kcheby,iv,j,ielec,iom)=0.d0
  else
- prodCvj(kcheby,iv,j,ielec,iom)=Cvjprod(iv,j,iom)
+ prodCvj(kcheby,iv,j,ielec,iom)=Cvjprod(iv,j,ican)
  endif
  enddo
  enddo

SRC/main_madwave3.f

@@ -527,15 +527,18 @@ subroutine printprob(iit0,iloop,indt)
  implicit none
  include "mpif.h"
  integer :: ifile,iii,ielec,ivprod,iomprod,ican,iele,j,iv
- integer :: ie,iom
+ integer :: ie,iom,ifilelec,ivfile,iiielec,iiiv
  integer :: iit0,iloop,indt
- real*8 :: S2reac,Av,S2no
+ real*8 :: S2reac,Av,S2no,S2nobis
  real*8 :: vibprod(nviniprod:nvmaxprod)
  real*8 :: rotdistri(jini:jmax)
 
  real*8 :: S2pro(nviniprod:nvmaxprod,jiniprod:jmaxprod
  & ,iomminprod:iommaxprod)
  complex*16 :: zzz
+
+ ifilelec=50
+ ivfile=100
 
 * writes information at each loop in time
 
@@ -544,10 +547,27 @@ subroutine printprob(iit0,iloop,indt)
  & nvref,Jtot,iloop
  open(20,file=name,status='unknown')
  else
- write(name,'("S2mat.J",i3.3,".k",i5.5)')
- & Jtot,iloop
- open(20,file=name,status='unknown')
+! write(name,'("S2mat.J",i3.3,".k",i5.5)')
+! & Jtot,iloop
+! open(20,file=name,status='unknown')
+
+ do ielec=1,nelec
+ iiielec=ifilelec+ielec
+ write(name,"('distriS2reac.elec',i2.2)")ielec
+ open(iiielec,file=name,status='unknown')
+ enddo
+ iiiv=ivfile
+ do ielec=1,nelec
+ do iv=nvini,noBCstates(jini,ielec)
+ iiiv=ivfile+(iv+1-nvini)+(nvmax-nvini+1)*(ielec-1)
+ write(name,"('distriS2reac.v',i2.2,'.e',i1.1)")
+ & iv,ielec
+ open(iiiv,file=name,status='unknown')
+ enddo
+ enddo
+
  endif
+
  do ie=1,netot
 
 ****> state-2-state for reactants
@@ -626,14 +646,53 @@ subroutine printprob(iit0,iloop,indt)
  endif
  else ! npun1 rotdis
 
- write(20,"(501(1x,e15.7))")etotS2(ie)/conve1/8065.5d0
- & ,S2reac*photonorm
+! write(20,"(501(1x,e15.7))")etotS2(ie)/conve1/8065.5d0
+! & ,S2reac*photonorm
+! & ,(rotdistri(j)*photonorm,j=jini,jmax)
+
+ iiiv=ivfile
+ do ielec=1,nelec
+ S2nobis=0.d0
+ do iv=nvini,noBCstates(jini,iele)
+ iiiv=ivfile+(iv+1-nvini)+(nvmax-nvini+1)*(ielec-1)
+ rotdistri(:)=0.d0
+ do ican=1,ncan
+ iele=nelebas(ican)
+ if(iele.eq.ielec)then
+ do j=j00(ican),jmax,inc
+ zzz=zS2(ie,iv,j,ican)/2.d0/pi
+ Av=dreal(zzz*dconjg(zzz))
+ rotdistri(j)=rotdistri(j)+S2(iv,j,ican)
+ S2nobis=S2nobis+S2(iv,j,ican)
+ enddo ! j
+ endif ! iele=ielec
+ enddo ! ican
+ write(iiiv,"(1000(1x,e15.7))")etotS2(ie)/conve1/8065.5d0
  & ,(rotdistri(j)*photonorm,j=jini,jmax)
 
+ enddo ! iv
+ iiielec=ifilelec+ielec
+ write(iiielec,"(1000(1x,e15.7))")
+ & etotS2(ie)/conve1/8065.5d0 ,S2nobis*photonorm
+
+ enddo ! ielec
+
  endif
  enddo ! ie=1,ne
-
- close(20)
+
+! close(20)
+
+ do ielec=1,nelec
+ iiielec=ifilelec+ielec
+ close(iiielec)
+ enddo
+ iiiv=ivfile
+ do ielec=1,nelec
+ do iv=nvini,noBCstates(jini,ielec)
+ iiiv=ivfile+(iv+1-nvini)+(nvmax-nvini+1)*(ielec-1)
+ close(iiiv)
+ enddo
+ enddo
 
  return
  end