diff --git a/src/dispersion.f90 b/src/dispersion.f90 index 85defcc..7cec6d0 100644 --- a/src/dispersion.f90 +++ b/src/dispersion.f90 @@ -3,9 +3,9 @@ module dispersion use calcei_mod implicit none integer, parameter :: dp=kind(1.0d0) - integer, parameter :: tmax=5,npts=500 + integer, parameter :: npts=500 real(dp), save :: ttv(npts+1),extv(npts+1) - real(dp), parameter :: dtex=2.0d0*tmax/dble(npts) + real(dp), parameter :: tmax=5.0d0,dtex=2.0d0*tmax/dble(npts) ! contains ! @@ -21,16 +21,14 @@ subroutine colddisp(xg,yg,npl,nprf,sox) real(dp) :: nprf ! N perpendicular to B (cold) real(dp) :: sox ! sox=-1 ==> O mode, sox=1 ==> X mode ! - real(dp) :: yg2,npl2,nprf2,n_f2,dnl,del + real(dp) :: yg2,npl2,dnl,del,dxg ! npl2=npl*npl dnl=1.0d0-npl2 + dxg=1.0d0-xg yg2=yg*yg - del=sqrt(dnl*dnl+4.0d0*npl2*(1.0d0-xg)/yg2) - n_f2=1.0d0-xg-xg*yg2*(1.0d0+npl2+sox*del)/(1.0d0-xg-yg2)/2.0d0 -! - nprf2=n_f2-npl2 - nprf=sqrt(nprf2) + del=sqrt(dnl*dnl+4.0d0*npl2*dxg/yg2) + nprf=sqrt(dxg-npl2-xg*yg2*(1.0d0+npl2+sox*del)/(dxg-yg2)/2.0d0) ! end subroutine colddisp ! @@ -403,7 +401,6 @@ subroutine diel_tens_fr(xg,yg,mu,npl,e330,cr,ci,epsl,lrm,fast) ! llm=min(3,lrm) ! - dt=2.0d0*tmax/dble(npts) bth2=2.0d0/mu bth=sqrt(bth2) mu2=mu*mu @@ -417,7 +414,7 @@ subroutine diel_tens_fr(xg,yg,mu,npl,e330,cr,ci,epsl,lrm,fast) upl2=bth2*x**2 upl=bth*x gx=1.0d0+t*t/mu - exdx=cr*extv(i)*gx/rxt*dt + exdx=cr*extv(i)*gx/rxt*dtex ! n1=1 if(fast.gt.2) n1=-llm