advection.f90 1.61 KB
!-------------------------------------------------------------------
!DESCRIPTION: Here, the wave spectrum is advected through the
!             domain. The advection equation is solved using a
!             Lax-wendroff discretization sheme with a Superbee
!             flux limiter.
!-------------------------------------------------------------------

!INTERFACE:
subroutine advection

   !MODULE USES:
   use parameters

   !LOCAL PARAMETERS:
   implicit none
                
   double precision, allocatable :: diffF(:)
   double precision, allocatable :: diffl(:)
   double precision, allocatable :: diffr(:)
   double precision, allocatable :: phi  (:)
   double precision, allocatable :: theta(:)
   double precision, allocatable :: F    (:)
   double precision, allocatable :: diff1(:)
        
   allocate(diffl(nx  ))
   allocate(diffr(nx  ))
   allocate(diffF(nx  ))
   allocate(theta(nx  ))
   allocate(phi  (nx  ))
   allocate(F    (nx  ))
   allocate(diff1(nx-1))

   if ( cont.eq.1 ) then
      E(1,ii,n-1)=Ei(ii)
   end if

   do i=2,nx
      diffl(i) = E(i  ,ii,n-1) - E(i-1,ii,n-1)  
   end do

   do i=1,nx-1
      diffr(i) = E(i+1,ii,n-1) - E(i  ,ii,n-1)
   end do

   ! set dE=0 at boundaries
   diffr(nx)=-E(nx,ii,n-1)
   diffl(1 )= E( 1,ii,n-1)

   ! Superbee flux limiter
   do i=1,nx
      theta(i)=diffl(i)/(diffr(i)+3e-14)
      phi(i)=max(0d0,min(theta(i)*2,1d0),min(theta(i),2d0))
   end do
             
   ! Lax-Wendroff sheme
   F=E(1:nx,ii,n-1) + 0.5*(1-CN(ii))*diffr*phi

   do i=2,nx
      diffF(i)=F(i) - F(i-1)
   end do

   diffF(1)=F(1)
   E(1:nx,ii,n)=E(1:nx,ii,n-1)-CN(ii)*diffF
   
end subroutine advection