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
- Jan S. Hesthaven and Tim Warburton. Nodal Discontinuous Galerkin Methods: Algorithms, Analysis, and Applications, Springer, 2008.
- The third-order Runge Kutta Method. Link.
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