Transparent boundary condition and its high frequency approximation for the Schrödinger equation on a rectangular computational domain

TL;DR

This work tackles nonreflecting boundary treatment for the 2D Schrödinger equation on a rectangle by comparing a high-frequency local boundary map with a novel Padé-based TBC. It develops a Legendre-Galerkin spectral discretization using Lobatto polynomials and couples it to time discretizations via convolution quadrature and Padé approximants to handle nonlocal temporal operators. The main contributions are a comprehensive HF boundary scheme with corner conditions, a storage-efficient Padé approach for both edge and corner TBCs, and extensive numerical experiments confirming stability and empirical convergence. The methods offer scalable, higher-order time stepping and enable rigorous comparison between HF and Padé-based TBCs, with potential applications to more general exterior potentials and nonlinear Schrödinger dynamics.

Abstract

This paper addresses the numerical implementation of the transparent boundary condition (TBC) and its various approximations for the free Schrödinger equation on a rectangular computational domain. In particular, we consider the exact TBC and its spatially local approximation under high frequency assumption along with an appropriate corner condition. For the spatial discretization, we use a Legendre-Galerkin spectral method where Lobatto polynomials serve as the basis. Within variational formalism, we first arrive at the time-continuous dynamical system using spatially discrete form of the initial boundary-value problem incorporating the boundary conditions. This dynamical system is then discretized using various time-stepping methods, namely, the backward-differentiation formula of order 1 and 2 (i.e., BDF1 and BDF2, respectively) and the trapezoidal rule (TR) to obtain a fully discrete system. Next, we extend this approach to the novel Padé based implementation of the TBC presented by Yadav and Vaibhav [arXiv:2403.07787(2024)]. Finally, several numerical tests are presented to demonstrate the effectiveness of the boundary maps (incorporating the corner conditions) where we study the stability and convergence behavior empirically.

Figures (8)

• Figure 1: The figure shows a rectangular domain with boundary segments parallel to one of the axes.
• Figure 2: A schematic depiction of the evolution of the auxiliary field $\varphi(x_1,x_2,\tau_1,\tau_2)$ in the $(\tau_1,\tau_2)$-plane is provided in this figure where the plot on the right corresponds $\boldsymbol{x}\in\Gamma_r\cup\Gamma_l$ and the plot on the left corresponds $\boldsymbol{x}\in\Gamma_t\cup\Gamma_b$. The filled circles depict the evolution of the auxiliary field $\varphi(x_1,x_2,\tau_1,\tau_2)$ either above or below the diagonal in the $(\tau_1,\tau_2)$-plane starting from the diagonal which also serves as initial conditions for solving IVPs corresponding to the auxiliary function. The TBCs for the auxiliary field require the history of the auxiliary field at the corner points which makes empty circles relevant. Note that these values at the corners can be taken from the adjacent segment of the boundary where it is already being computed and this is depicted by broken lines. Note that the vertical/horizontal lines where the arrows end corresponds to the history of the auxiliary field needed for the TBCs on $\partial\Omega_i$ in the current time ($t$).
• Figure 3: A schematic depiction of the evolution of the auxiliary fields $\varphi_{k,a_1}(x_2,\tau_1,\tau_2),\; \varphi_{k,a_2}(x_1,\tau_1,\tau_2)$ and $\psi_{k,k',a_1,a_2}(\tau_1,\tau_2)$ in the $(\tau_1,\tau_2)$-plane is provided in this figure. The plots (A) and (B) depict the evolution of the fields $\varphi_{k,a_2}(x_1,\tau_1,\tau_2)$ and $\varphi_{k,a_1}(x_2,\tau_1,\tau_2)$ on the boundary segments $\Gamma_{a_2}$ and $\Gamma_{a_1}$, respectively. The plot (C) depicts the evolution of the field $\psi_{k,k',a_1,a_2}(\tau_1,\tau_2)$ which can be achieved by moving either below or above the diagonal.
• Figure 4: The figure shows the evolution of the relative energy content as defined in \ref{['eq:cg2d-energy-content']} of the chirped-Gaussian and Hermite-Gaussian profiles considered in Table \ref{['tab:cg2d']} and Table \ref{['tab:hg2d']}. Here the computational domain is $\Omega_i=(-10,10)^2$.
• Figure 5: The figure shows a comparison of evolution of error in the numerical solution of the IBVP \ref{['eq:2D-SE-CT']} with various approximations of the TBCs for the chirped-Gaussian and Hermite-Gaussian profiles with different values of the speed 'c' (see Table \ref{['tab:cg2d']} and Table \ref{['tab:hg2d']}). The numerical parameters and the labels are described in Sec. \ref{['sec:tests-ee']} where the error is quantified by \ref{['eq:error-ibvp']}.

Theorems & Definitions (4)

• Remark 1
• Remark 2
• Remark 3
• Remark 4