Sunday, June 24, 2012

LAKE: DG-Fourier Multi-stage SSP method for the 2-D advection equation


Purpose
  • The package "Lake" solves the 2-D advection equation in a rectangular domain. A Fourier psuedospectral method was used for the spatial discretization in $x$, and a nodal discontinuous Galerkin method was used for $y$. The temporal discretization are one-step multi-stage explicit schemes. It is partially based on packages "Sea" and "Ion".

Specifications
  • Name: Lake.
  • Author: Feng Chen.
  • Finishing date: 06/26/2012.
  • Languages: Fortran 90.
  • Required libraries: BLAS, LAPACK.

Simple Example
  • Equation: \begin{equation} \begin{aligned} & \Psi_t + \boldsymbol{u} \cdot \nabla \Psi = f, \\ & \Psi(x,y,0) = \Psi_0(x,y), \\ & \Psi \text{ periodic in } x, \\ & \Psi \text{ no boundary condition in } y, \end{aligned} \end{equation} where $\boldsymbol{u}$ satisfies \begin{equation} \nabla \cdot \boldsymbol{u} = 0, \quad \boldsymbol{u} \cdot \boldsymbol{n} |_{y=\pm 1}= 0. \end{equation}
  • Parameters: \begin{equation} \Omega =(-\pi,\pi)\times (-1,1), \, T=10, \, \delta t = 0.001, \, N_x = 64, \, N_y = 10, \, K_y = 10. \end{equation}
  • Exact solution and input functions: \begin{equation} \begin{aligned} & \Psi(x,y,t) = \cos(x+t) \sin(y+t),\\ & u_1(x,y,t) = \pi \cos(x+t) \cos(\pi y), \\ &u_2(x,y,t) = \sin(x+t) \sin(\pi y), \end{aligned} \end{equation} then $f(x,y,t)$ is calculated accordingly.

Quick Start
  • Compiling and running:
    cd ./Lake
    make library
    gfortran Lake_Main.cu -llibrary -llapack -lblas
    ./a.out
    
  • Output:
    Order 3 Runge Kutta was used.
     at step           0 max error=   0.0000000000000000     
     at step        2000 max error=  2.32462854521386930E-009
     at step        4000 max error=  2.56192570491364791E-008
     at step        6000 max error=  6.41112042915753522E-008
     at step        8000 max error=  1.24243067334273150E-007
     at step       10000 max error=  1.15325957505962862E-007
    ...
    
  • CPU: Intel(R) Xeon(R) CPU X5550 @2.67GHz.
  • OS: CentOS release 6.4 (Final).
  • Compiler: gfortran 4.5.1.

References
Code Highlights
    
    !<
    !! calculate the rhs of the linear system
    !<
    subroutine LakeSpaceTime_RHS(SI, SF, rhs, flux, Psi, fsrc)  
     type(domain), intent(in) ::  SF   
     type(IonSpace1D), intent(in) :: SI 
     real(kind=8), dimension(SI%n,SF%ndof), intent(in) :: Psi, fsrc
     real(kind=8), dimension(SI%n,SF%ndof) :: rhs
     real(kind=8), dimension(SI%n,SF%ndof,2) :: flux 
     integer :: j, k  
    
     ! differentiation in y
     do k = 1, SI%n
      call NDG_RHS(SF, rhs(k,:), flux(k,:,2), Psi(k,:))
     end do
     flux(:,:,2) = rhs
    
     ! differentiation in x
     do j = 1, SF%ndof
      ! transform to k-space
      call rfft1d(SI%n, flux(:,j,1), SI%wk, 0)
      ! differentiation
      call fodxcf(SI%n, flux(:,j,1), SI%sc)  
      ! transform to real space
      call rfft1d(SI%n, flux(:,j,1), SI%wk, 1)
     end do 
     ! pseudo-Fourier, no integration by parts.
     flux(:,:,1) = -flux(:,:,1)  
    
    
     ! add up all the rhs
       rhs = flux(:,:,1) + flux(:,:,2) + fsrc 
    end subroutine LakeSpaceTime_RHS
    
    

No comments:

Post a Comment