Craig C. Douglas, Long Lee, Man-Chung Yeung
This paper presents the first results to combine two theoretically sound methods (spectral projection and multigrid methods) together to attack ill-conditioned linear systems. Our preliminary results show that the proposed algorithm applied to a Krylov subspace method takes much fewer iterations for solving an ill-conditioned problem downloaded from a popular online sparse matrix collection.
Yi-Hung Kuo, Long Lee, Gregory Lyng
Finite difference approximation, in addition to Taylor truncation errors, introduces numerical dispersion-and-dissipation errors into numerical solutions of partial differential equations. We analyze a class of finite difference schemes which are designed to minimize these errors (at the expense of formal order of accuracy), and we analyze the interplay between the Taylor truncation errors and the dispersion-and-dissipation errors during mesh refinement. In particular, we study the numerical dispersion relation of the fully discretized non-dispersive transport equation in one and two space dimensions. We derive the numerical phase error and the $L^2$-norm error of the solution in terms of the dispersion-and-dissipation error. Based on our analysis, we investigate the error dynamics among various optimized compact schemes and the unoptimized higher-order generalized Padé compact schemes. The dynamics shed light on the principles of designing suitable optimized compact schemes for a given problem. Using these principles as guidelines, we then propose an optimized scheme that prescribes the numerical dispersion relation before finding the corresponding discretization. This approach produces smaller numerical dispersion-and-dissipation errors for linear and nonlinear problems, compared with the unoptimized higher-order compact schemes and other optimized schemes developed in the literature. Finally, we discuss the difficulty of developing an optimized composite boundary scheme for problems with non-trivial boundary conditions. We propose a composite scheme that introduces a buffer zone to connect an optimized interior scheme and an unoptimized boundary scheme. Our numerical experiments show that this strategy produces small $L^2$-norm error when a wave packet passes through the non-periodic boundary.
Gastão A. Braga, Federico C. Furtado, Vincenzo Isaia et al.
We systematically study a numerical procedure that reveals the asymptotically self-similar dynamics of solutions of partial differential equations (PDEs). This procedure, based on the renormalization group (RG) theory for PDEs, appeared initially in a conference proceeding by Braga et al. \cite{BFI04}. This numerical version of RG method, dubbed as the numerical RG (nRG) algorithm, numerically rescales the temporal and spatial variables in each iteration and drives the solutions to a fixed point exponentially fast, which corresponds to the self-similar dynamics of the equations. In this paper, we carefully examine and validate this class of algorithms by comparing the numerical solutions with either the exact or the asymptotic solutions of the model equations in literature. The other contribution of the current paper is that we present several examples to demonstrate that this class of nRG algorithms can be applied to a wide range of PDEs to shed lights on longtime self-similar dynamics of certain physical models that are difficult to analyze, both numerically and analytically.
Man-Chung Yeung, Craig C. Douglas, Long Lee
We present a preconditioner based on spectral projection that is combined with a deflated Krylov subspace method for solving ill conditioned linear systems of equations. Our results show that the proposed algorithm requires many fewer iterations to achieve the convergence criterion for solving an ill conditioned problem than a Krylov subspace solver. In our numerical experiments, the solution obtained by the proposed algorithm is more accurate in terms of the norm of the distance to the exact solution of the linear system of equations.
S. Huzurbazar, Long Lee, Dongyang Kuang
We study a class of mathematical and statistical algorithms with the aim of establishing a computer-based framework for fast and reliable automatic abnormality detection on landmark represented image templates. Under this framework, we apply a landmark-based algorithm for finding a group average as an estimator that is said to best represent the common features of the group in study. This algorithm extracts information of momentum at each landmark through the process of template matching. If ever converges, the proposed algorithm produces a local coordinate system for each member of the observing group, in terms of the residual momentum. We use a Bayesian approach on the collected residual momentum representations for making inference. For illustration, we apply this framework to a small database of brain images for detecting structure abnormality. The brain structure changes identified by our framework are highly consistent with studies in the literature.
Roberto Camassa, Dongyang Kuang, Long Lee
The Euler-Poincaré (EP) equations describe the geodesic motion on the diffeomorphism group. For template matching (template deformation), the Euler-Lagrangian equation, arising from minimizing an energy function, falls into the Euler-Poincaré theory and can be recast into the EP equations. By casting the EP equations in the Lagrangian (or characteristics) form, we formulate the equations as a finite dimensional particle system. The evolution of this particle system describes the geodesic motion of landmark points on a Riemann manifold. In this paper we present a class of novel algorithms that take advantage of the structure of the particle system to achieve a fast matching process between the reference and the target templates. The strong suit of the proposed algorithms includes (1) the efficient feedback control iteration, which allows one to find the initial velocity field for driving the deformation from the reference template to the target one, (2) the use of the conical kernel in the particle system, which limits the interaction between particles and thus accelerates the convergence, and (3) the availability of the implementation of fast-multipole method for solving the particle system, which could reduce the computational cost from $O(N^2)$ to $O(N\log N)$, where $N$ is the number of particles. The convergence properties of the proposed algorithms are analyzed. Finally, we present several examples for both exact and inexact matchings, and numerically analyze the iterative process to illustrate the efficiency and the robustness of the proposed algorithms.
Roberto Camassa, Dongyang Kuang, Long Lee
We study a class of partial differential equations (PDEs) in the family of the so-called Euler-Poincaré differential systems, with the aim of developing a foundation for numerical algorithms of their solutions. This requires particular attention to the mathematical properties of this system when the associated class of elliptic operators possesses non-smooth kernels. By casting the system in its Lagrangian (or characteristics) form, we first formulate a particles system algorithm in free space with homogeneous Dirichlet boundary conditions for the evolving fields. We next examine the deformation of the system when non-homogeneous "constant stream" boundary conditions are assumed. We show how this simple change at the boundary deeply affects the nature of the evolution, from hyperbolic-like to dispersive with a non-trivial dispersion relation, and examine the potentially regularizing properties of singular kernels offered by this deformation. From the particle algorithm viewpoint, kernel singularities affect the existence and uniqueness of solutions to the corresponding ordinary differential equations systems. We illustrate this with the case when the operator kernel assumes a conical shape over the spatial variables, and examine in detail two-particle dynamics under the resulting lack of Lipschitz-continuity. Curiously, we find that for the conically-shaped kernels the motion of the related two-dimensional waves can become completely integrable under appropriate initial data. This reduction projects the two-dimensional system to the one-dimensional completely integrable Shallow-Water equation [Camassa, R. and Holm, D. D., Phys. Rev. Lett., 71, 1961-1964, 1993], while retaining the full dependence on two spatial dimensions for the single channel solutions.
Dongyang Kuang, Long Lee
We present a conservation formulation and a numerical algorithm for the reduced-gravity shallow-water equations on a beta plane, subjected to a constant wind forcing that leads to the formation of double-gyre circulation in a closed ocean basin. The novelty of the paper is that we reformulate the governing equations into a nonlinear hyperbolic conservation law plus source terms. A second-order fractional-step algorithm is used to solve the reformulated equations. In the first step of the fractional-step algorithm, we solve the homogeneous hyperbolic shallow-water equations by the wave-propagation finite volume method. The resulting intermediate solution is then used as the initial condition for the initial-boundary value problem in the second step. As a result, the proposed method is not sensitive to the choice of viscosity and gives high-resolution results for coarse grids, as long as the Rossby deformation radius is resolved. We discuss the boundary conditions in each step, when no-slip boundary conditions are imposed to the problem. We validate the algorithm by a periodic flow on an f-plane with exact solutions. The order-of-accuracy for the proposed algorithm is tested numerically. We illustrate a quasi-steady-state solution of the double-gyre model via the height anomaly and the contour of stream function for the formation of double-gyre circulation in a closed basin. Our calculations are highly consistent with the results reported in the literature. Finally, we present an application, in which the double-gyre model is coupled with the advection equation for modeling transport of a pollutant in a closed ocean basin.