Computation of the inverse Laplace transform based on a collocation method which uses only real values

We develop a numerical algorithm for inverting a Laplace transform (LT), based on Laguerre polynomial series expansion of the inverse function under the assumption that the LT is known on the real axis only. The method belongs to the class of Collocation methods (C-methods), and is applicable when the LT function is regular at infinity. Difficulties associated with these problems are due to their intrinsic ill-posedness. The main contribution of this paper is to provide computable estimates of truncation, discretization, conditioning and roundoff errors introduced by numerical computations. Moreover, we introduce the pseudoaccuracy which will be used by the numerical algorithm in order to provide uniform scaled accuracy of the computed approximation for any x   with respect to e^σx${\mathrm{e}}^{\sigma x}$. These estimates are then employed to dynamically truncate the series expansion. In other words, the number of the terms of the series acts like the regularization parameter which provides the trade-off between errors.     

On the condition number of matrices arising in the tikhonov regularization method

In this work, we consider the regularization method for linear ill-posed problems. For operators and approximating subspaces satisfying certain conditions and for a specific form of the regularization parameter, upper and lower bounds are given for the condition number of the corresponding discrete problem.     

On deflation and singular symmetric positive semi-definite matrices

For various applications, it is well-known that the deflated ICCG is an efficient method for solving linear systems with invertible coefficient matrix. We propose two equivalent variants of this deflated ICCG which can also solve linear systems with singular coefficient matrix, arising from discretization of the discontinuous Poisson equation with Neumann boundary conditions. It is demonstrated both theoretically and numerically that the resulting methods accelerate the convergence of the iterative process.     

On algebraic multi-level methods for non-symmetric systems - Comparison results

We establish theoretical comparison results for algebraic multi-level methods applied to non-singular non-symmetric M-matrices. We consider two types of multi-level approximate block factorizations or AMG methods, the AMLI and the MAMLI method. We compare the spectral radii of the iteration matrices of these methods. This comparison shows, that the spectral radius of the MAMLI method is less than or equal to the spectral radius of the AMLI method. Moreover, we establish how the quality of the approximations in the block factorization effects the spectral radii of the iteration matrices. We prove comparisons results for different approximations of the fine grid block as well as for the used Schur complement. We also establish a theoretical comparison between the AMG methods and the classical block Jacobi and block Gauss-Seidel methods.     

On nonsymmetric saddle point matrices that allow conjugate gradient iterations

Linear systems in saddle point form are usually highly indefinite,which often slows down iterative solvers such as Krylov subspace methods. It has been noted by several authors that negating the second block row of a symmetric indefinite saddle point matrix leads to a nonsymmetric matrix ${{\mathcal A}}Linear systems in saddle point form are usually highly indefinite,which often slows down iterative solvers such as Krylov subspace methods. It has been noted by several authors that negating the second block row of a symmetric indefinite saddle point matrix leads to a nonsymmetric matrix whose spectrum is entirely contained in the right half plane. In this paper we study conditions so that is diagonalizable with a real and positive spectrum. These conditions are based on necessary and sufficient conditions for positive definiteness of a certain bilinear form,with respect to which is symmetric. In case the latter conditions are satisfied, there exists a well defined conjugate gradient (CG) method for solving linear systems with . We give an efficient implementation of this method, discuss practical issues such as error bounds, and present numerical experiments. In memory of Gene Golub (1932–2007), our wonderful friend and colleague, who had a great interest in the conjugate gradient method and the numerical solution of saddle point problems. The work of J?rg Liesen was supported by the Emmy Noether-Program and the Heisenberg-Program of the Deutsche Forschungsgemeinschaft.     

On CSCS-based iteration methods for Toeplitz system of weakly nonlinear equations

For Toeplitz system of weakly nonlinear equations, by using the separability and strong dominance between the linear and the nonlinear terms and using the circulant and skew-circulant splitting (CSCS) iteration technique, we establish two nonlinear composite iteration schemes, called Picard-CSCS and nonlinear CSCS-like iteration methods, respectively. The advantage of these methods is that they do not require accurate computation and storage of Jacobian matrix, and only need to solve linear sub-systems of constant coefficient matrices. Therefore, computational workloads and computer storage may be saved in actual implementations. Theoretical analysis shows that these new iteration methods are local convergent under suitable conditions. Numerical results show that both Picard-CSCS and nonlinear CSCS-like iteration methods are feasible and effective for some cases.     

The truncated Stieltjes moment problem solved by using kernel density functions

In this work the problem of the approximate numerical determination of a semi-infinite supported, continuous probability density function (pdf) from a finite number of its moments is addressed. The target space is carefully defined and an approximation theorem is proved, establishing that the set of all convex superpositions of appropriate Kernel Density Functions (KDFs) is dense in this space. A solution algorithm is provided, based on the established approximate representation of the target pdf and the exploitation of some theoretical results concerning moment sequence asymptotics. The solution algorithm also permits us to recover the tail behavior of the target pdf and incorporate this information in our solution. A parsimonious formulation of the proposed solution procedure, based on a novel sequentially adaptive scheme is developed, enabling a very efficient moment data inversion. The whole methodology is fully illustrated by numerical examples.     

On the regularizing properties of the GMRES method

Summary. The GMRES method is a popular iterative method for the solution of large linear systems of equations with a nonsymmetric nonsingular matrix. However, little is known about the behavior of this method when it is applied to the solution of nonsymmetric linear ill-posed problems with a right-hand side that is contaminated by errors. We show that when the associated error-free right-hand side lies in a finite-dimensional Krylov subspace, the GMRES method is a regularization method. The iterations are terminated by a stopping rule based on the discrepancy principle. Received November 10, 2000 / Revised version received April 11, 2001 / Published online October 17, 2001     

Preconditioning of elliptic problems by approximation in the transform domain

Preconditioned conjugate gradient method is applied for solving linear systemsAx=b where the matrixA is the discretization matrix of second-order elliptic operators. In this paper, we consider the construction of the trnasform based preconditioner from the viewpoint of image compression. Given a smooth image, a major portion of the energy is concentrated in the low frequency regions after image transformation. We can view the matrixA as an image and construct the transform based preconditioner by using the low frequency components of the transformed matrix. It is our hope that the smooth coefficients of the given elliptic operator can be approximated well by the low-rank matrix. Numerical results are reported to show the effectiveness of the preconditioning strategy. Some theoretical results about the properties of our proposed preconditioners and the condition number of the preconditioned matrices are discussed.     

On the eigenvalues of a class of saddle point matrices

We study spectral properties of a class of block 2 × 2 matrices that arise in the solution of saddle point problems. These matrices are obtained by a sign change in the second block equation of the symmetric saddle point linear system. We give conditions for having a (positive) real spectrum and for ensuring diagonalizability of the matrix. In particular, we show that these properties hold for the discrete Stokes operator, and we discuss the implications of our characterization for augmented Lagrangian formulations, for Krylov subspace solvers and for certain types of preconditioners. The work of this author was supported in part by the National Science Foundation grant DMS-0207599 Revision dated 5 December 2005.     

Numerical accuracy of real inversion formulas for the Laplace transform

In this article, we investigate and compare a number of real inversion formulas for the Laplace transform. The focus is on the accuracy and applicability of the formulas for numerical inversion. In this contribution, we study the performance of the formulas for measures concentrated on a positive half-line to continue with measures on an arbitrary half-line. As our trial measure concentrated on a positive half-line, we take the broad Gamma probability distribution family.     

A preconditioning approach to the pagerank computation problem

Some spectral properties of the transition matrix of an oriented graph indicate the preconditioning of Euler-Richardson (ER) iterative scheme as a good way to compute efficiently the vertexrank vector associated with such graph. We choose the preconditioner from an algebra U of matrices, thereby obtaining an ER_U method, and we observe that ER_U can outperform ER in terms of rate of convergence. The proposed preconditioner can be updated at a very low cost whenever the graph changes, as is the case when it represents a generic set of information. The particular U utilized requires a surplus of operations per step and memory allocations, which make ER_U superior to ER for not too wide graphs. However, the observed high improvement in convergence rate obtained by preconditioning and the general theory developed, are a reason for investigating different choices of U, more efficient for huge graphs.     

Computing the Hilbert transform of a Jacobi weight function

We discuss the evaluation of the Hilbert transformf _–1 ¹ (t-)^–1 w(, )(t)dt,–1<<1, of the Jacobi weight functionw(, )(t)=(1–t))(1+t) by analytic and numerical means and also comment on the recursive computation of the quantitiesf _–1 ¹ )(t–)^–1 _n (t;w ^{(, ))} w ^{(, )}(t)dt,n=0, 1, 2, ..., where _n (·;w ^{(, )}) is the Jacobi polynomial of degreen.The work of the first author was supported in part by the National Science Foundation under grant DCR-8320561. The work of the second author was supported by the National Science Foundation under grant DMS-8419086.     

On the modification of an eigenvalue problem that preserves an eigenspace

Eigenvalue problems arise in many application areas ranging from computational fluid dynamics to information retrieval. In these fields we are often interested in only a few eigenvalues and corresponding eigenvectors of a sparse matrix. In this paper, we comment on the modifications of the eigenvalue problem that can simplify the computation of those eigenpairs. These transformations allow us to avoid difficulties associated with non-Hermitian eigenvalue problems, such as the lack of reliable non-Hermitian eigenvalue solvers, by mapping them into generalized Hermitian eigenvalue problems. Also, they allow us to expose and explore parallelism. They require knowledge of a selected eigenvalue and preserve its eigenspace. The positive definiteness of the Hermitian part is inherited by the matrices in the generalized Hermitian eigenvalue problem. The position of the selected eigenspace in the ordering of the eigenvalues is also preserved under certain conditions. The effect of using approximate eigenvalues in the transformation is analyzed and numerical experiments are presented.     

Subspace-restricted singular value decompositions for linear discrete ill-posed problems

The truncated singular value decomposition is a popular solution method for linear discrete ill-posed problems. These problems are numerically underdetermined. Therefore, it can be beneficial to incorporate information about the desired solution into the solution process. This paper describes a modification of the singular value decomposition that permits a specified linear subspace to be contained in the solution subspace for all truncations. Modifications that allow the range to contain a specified subspace, or that allow both the solution subspace and the range to contain specified subspaces also are described.     

Multi-dimensional option pricing using radial basis functions and the generalized Fourier transform

We show that the generalized Fourier transform can be used for reducing the computational cost and memory requirements of radial basis function methods for multi-dimensional option pricing. We derive a general algorithm, including a transformation of the Black–Scholes equation into the heat equation, that can be used in any number of dimensions. Numerical experiments in two and three dimensions show that the gain is substantial even for small problem sizes. Furthermore, the gain increases with the number of dimensions.     

Vector coding algorithms for multidimensional discrete Fourier transform

A new fast algorithm is presented for the multidimensional discrete Fourier transform (DFT). This algorithm is derived using an interesting technique called “vector coding” (VC), and we call it the vector-coding fast Fourier transform (VC-FFT) algorithm. Since the VC-FFT is an extension of the Cooley–Tukey algorithm from 1-D to multidimensional form, the structure of the program is as simple as the Cooley–Tukey fast Fourier transform (FFT). The new algorithm significantly reduces the number of multiplications and recursive stages. The VC-FFT therefore comprehensively reduces the complexity of the algorithm as compared with other current multidimensional DFT algorithms.     

Estimating the largest singular values of large sparse matrices via modified moments

We describe a procedure for determining a few of the largest singular values of a large sparse matrix. The method by Golub and Kent which uses the method of modified moments for estimating the eigenvalues of operators used in iterative methods for the solution of linear systems of equations is appropriately modified in order to generate a sequence of bidiagonal matrices whose singular values approximate those of the original sparse matrix. A simple Lanczos recursion is proposed for determining the corresponding left and right singular vectors. The potential asynchronous computation of the bidiagonal matrices using modified moments with the iterations of an adapted Chebyshev semi-iterative (CSI) method is an attractive feature for parallel computers. Comparisons in efficiency and accuracy with an appropriate Lanczos algorithm (with selective re-orthogonalization) are presented on large sparse (rectangular) matrices arising from applications such as information retrieval and seismic reflection tomography. This procedure is essentially motivated by the theory of moments and Gauss quadrature.This author's work was supported by the National Science Foundation under grants NSF CCR-8717492 and CCR-910000N (NCSA), the U.S. Department of Energy under grant DOE DE-FG02-85ER25001, and the Air Force Office of Scientific Research under grant AFOSR-90-0044 while at the University of Illinois at Urbana-Champaign Center for Supercomputing Research and Development.This author's work was supported by the U.S. Army Research Office under grant DAAL03-90-G-0105, and the National Science Foundation under grant NSF DCR-8412314.     

On generalized symmetric SOR method for augmented systems

In this paper, we establish the generalized symmetric SOR method (GSSOR) for solving the large sparse augmented systems of linear equations, which is the extension of the SSOR iteration method. The convergence of the GSSOR method for augmented systems is studied. Numerical resume shows that this method is effective.