Search

Transition matrices are a tool used to describe transition probabilities of a Markovian process. From an applied perspective, transition matrices are, for example, used in the fields of credit risk modeling and medical decision analysis. In these domains, it is possible to estimate matrices according to the frequency of data observations, but for modeling purposes, it may be necessary to calculate transition matrices for shorter and more restricted time intervals. From a mathematical perspective, this requirement translates into calculating the pth root of a stochastic matrix through a constrained optimization problem. This article focuses on computing pth root of transition matrices and proposes a novel approach using deep learning-based algorithms, particularly Convolutional Neural Networks. Comparisons with traditional algorithms reveal that, in certain contexts, deep learning models offer computational advantages, especially regarding time efficiency. The applicability of the models is verified through a case study built on Fitch Ratings data. [ABSTRACT FROM AUTHOR]

Algorithms for the computation of the (weighted) geometric mean G of two positive definite matrices are described and discussed. For large and sparse matrices the problem of computing the product y = G b , and of solving the linear system G x = b , without forming G, is addressed. An analysis of the conditioning is provided. Substantial numerical experimentation is carried out to test and compare the performances of these algorithms in terms of CPU time, numerical stability, and number of iterative steps. [ABSTRACT FROM AUTHOR]

This paper deals with the computation of the two parameter Mittag-Leffler function of operators by exploiting its Stieltjes integral representation and then by using a single exponential transform together with the sinc rule. Whenever the parameters of the function do not allow this representation, we resort to the Dunford-Taylor one. The error analysis is kept in the framework of unbounded accretive operators in order to make it a useful tool for the solution of fractional differential equations. The theory is also used to design a rational Krylov method. [ABSTRACT FROM AUTHOR]

This article considers the computation of the matrix exponential eA{{\rm{e}}}^{A} with numerical quadrature. Although several quadrature-based algorithms have been proposed, they focus on (near) Hermitian matrices. In order to deal with non-Hermitian matrices, we use another integral representation including an oscillatory term and consider applying the double exponential (DE) formula specialized to Fourier integrals. The DE formula transforms the given integral into another integral whose interval is infinite, and therefore, it is necessary to truncate the infinite interval. In this article, to utilize the DE formula, we analyze the truncation error and propose two algorithms. The first one approximates eA{{\rm{e}}}^{A} with the fixed mesh size, which is a parameter in the DE formula affecting the accuracy. The second one computes eA{{\rm{e}}}^{A} based on the first one with automatic selection of the mesh size depending on the given error tolerance.

It is well-known that the two-parameter Mittag-Leffler (ML) function plays a key role in Fractional Calculus. In this paper, we address the problem of computing this function, when its argument is a square matrix. Effective methods for solving this problem involve the computation of higher order derivatives or require the use of mixed precision arithmetic. In this paper, we provide an alternative method that is derivative-free and works entirely using IEEE standard double precision arithmetic. If certain conditions are satisfied, our method uses a Taylor series representation for the ML function; if not, it switches to a Schur-Parlett technique that will be combined with the Cauchy integral formula. A detailed discussion on the choice of a convenient contour is included. Theoretical and numerical issues regarding the performance of the proposed algorithm are discussed. A set of numerical experiments shows that our novel approach is competitive with the state-of-the-art method for IEEE double precision arithmetic, in terms of accuracy and CPU time. For matrices whose Schur decomposition has large blocks with clustered eigenvalues, our method far outperforms the other. Since our method does not require the efficient computation of higher order derivatives, it has the additional advantage of being easily extended to other matrix functions (e.g., special functions). [ABSTRACT FROM AUTHOR]

Nonsymmetric algebraic Riccati equations are a class of nonlinear matrix equations with wide applications. In this paper, numerical solution of nonsymmetric algebraic Riccati equations is studied. A generalized ALI iteration method is proposed to solve the equation, which is an extension of the ALI iteration method and the MALI iteration method. Theoretical analysis and numerical experiments show that the proposed method is feasible and is more effective than the ALI iteration method and the MALI iteration method. [ABSTRACT FROM AUTHOR]

To integrate large systems of nonlinear differential equations in time, we consider a variant of nonlinear waveform relaxation (also known as dynamic iteration or Picard–Lindelöf iteration), where at each iteration a linear inhomogeneous system of differential equations has to be solved. This is done by the exponential block Krylov subspace (EBK) method. Thus, we have an inner-outer iterative method, where iterative approximations are determined over a certain time interval, with no time stepping involved. This approach has recently been shown to be efficient as a time-parallel integrator within the PARAEXP framework. In this paper, convergence behavior of this method is assessed theoretically and practically. We examine efficiency of the method by testing it on nonlinear Burgers, Liouville–Bratu–Gelfand, and nonlinear heat conduction equations and comparing its performance with that of conventional time-stepping integrators. [ABSTRACT FROM AUTHOR]

The tensor t-function, a formalism that generalizes the well-known concept of matrix functions to third-order tensors, is introduced in Lund (Numer Linear Algebra Appl 27(3):e2288). In this work, we investigate properties of the Fréchet derivative of the tensor t-function and derive algorithms for its efficient numerical computation. Applications in condition number estimation and nuclear norm minimization are explored. Numerical experiments implemented by the t-Frechet toolbox hosted at https://gitlab.com/katlund/t-frechet illustrate properties of the t-function Fréchet derivative, as well as the efficiency and accuracy of the proposed algorithms. [ABSTRACT FROM AUTHOR]

A new iterative method for solving large scale symmetric nonlinear eigenvalue problems is presented. We firstly derive an infinite dimensional symmetric linearization of the nonlinear eigenvalue problem, then we apply the indefinite Lanczos method to this specific linearization, resulting in a short-term recurrence. We show how, under specific assumption on the starting vector, this method can be carried out in finite arithmetic and how the exploitation of the problem structure leads to improvements in terms of computation time. The eigenpair approximations are extracted with the nonlinear Rayleigh-Ritz procedure combined with a specific choice of the projection space. We illustrate how this extraction technique resolves the instability issues that may occur due to the loss of orthogonality in many standard Lanczos-type methods. [ABSTRACT FROM AUTHOR]

The use of fractional differential equations is a key tool in modeling non-local phenomena. Often, an efficient scheme for solving a linear system involving the discretization of a fractional operator is computing inverse fractional powers of the standard discretized Laplace operator. In this work, an exponential sum approximation for such fractional powers is derived. It is accurate over all positive real numbers larger than one, and allows to efficiently approximate the action of such operators on tensors stored in a variety of low-rank formats (CP, TT, Tucker). The results are relevant from a practical and theoretical perspective, as they predict the low-rank approximability of the solutions of these linear systems in low-rank tensor formats. [ABSTRACT FROM AUTHOR]

We consider the approximation of the inverse square root of regularly accretive operators in Hilbert spaces. The approximation is of rational type and comes from the use of the Gauss–Legendre rule applied to a special integral formulation of the fractional power. We derive sharp error estimates, based on the use of the numerical range, and provide some numerical experiments. For practical purposes, the finite-dimensional case is also considered. In this setting, the convergence is shown to be of exponential type. The method is also tested for the computation of a generic fractional power. [ABSTRACT FROM AUTHOR]

Various methods have been introduced for calculating the matrix sine and cosine. This paper concerns the Fibonacci-Horner decomposition of the matrix sine and cosine sin(tA), cos(tA) (A ∈ M (r; ℂ), t ∈ ℝ), which is derived from the combinatorial properties of the generalized Fibonacci sequences in the algebra of square matrices. Finally we derive some formulas for the trigonometric functions of f(A) where f is an homographic function. [ABSTRACT FROM AUTHOR]

We investigate the problem of approximating the matrix function f(A) by r(A), with f a Markov function, r a rational interpolant of f, and A a symmetric Toeplitz matrix. In a first step, we obtain a new upper bound for the relative interpolation error 1 − r/f on the spectral interval of A. By minimizing this upper bound over all interpolation points, we obtain a new, simple, and sharp a priori bound for the relative interpolation error. We then consider three different approaches of representing and computing the rational interpolant r. Theoretical and numerical evidence is given that any of these methods for a scalar argument allows to achieve high precision, even in the presence of finite precision arithmetic. We finally investigate the problem of efficiently evaluating r(A), where it turns out that the relative error for a matrix argument is only small if we use a partial fraction decomposition for r following Antoulas and Mayo. An important role is played by a new stopping criterion which ensures to automatically find the degree of r leading to a small error, even in presence of finite precision arithmetic. [ABSTRACT FROM AUTHOR]

In the stiff situation, we consider the long-time behavior of the relative error γ n in the numerical integration of a linear ordinary differential equation y ′ (t) = A y (t) , t ≥ 0 , where A is a normal matrix. The numerical solution is obtained by using at any step an approximation of the matrix exponential, e.g. a polynomial or a rational approximation. We study the long-time behavior of γ n by comparing it to the relative error γ n long in the numerical integration of the long-time solution, i.e. the projection of the solution on the eigenspace of the rightmost eigenvalues. The error γ n long grows linearly in time, it is small and it remains small in the long-time. We give a condition under which γ n ≈ γ n long , i.e. γ n γ n long ≈ 1 , in the long-time. When this condition does not hold, the ratio γ n γ n long is large for all time. These results describe the long-time behavior of the relative error γ n in the stiff situation. [ABSTRACT FROM AUTHOR]

Randomized trace estimation is a popular and well-studied technique that approximates the trace of a large-scale matrix B by computing the average of x T B x for many samples of a random vector X. Often, B is symmetric positive definite (SPD) but a number of applications give rise to indefinite B. Most notably, this is the case for log-determinant estimation, a task that features prominently in statistical learning, for instance in maximum likelihood estimation for Gaussian process regression. The analysis of randomized trace estimates, including tail bounds, has mostly focused on the SPD case. In this work, we derive new tail bounds for randomized trace estimates applied to indefinite B with Rademacher or Gaussian random vectors. These bounds significantly improve existing results for indefinite B, reducing the number of required samples by a factor n or even more, where n is the size of B. Even for an SPD matrix, our work improves an existing result by Roosta-Khorasani and Ascher (Found Comput Math, 15(5):1187–1212, 2015) for Rademacher vectors. This work also analyzes the combination of randomized trace estimates with the Lanczos method for approximating the trace of f(B). Particular attention is paid to the matrix logarithm, which is needed for log-determinant estimation. We improve and extend an existing result, to not only cover Rademacher but also Gaussian random vectors. [ABSTRACT FROM AUTHOR]

In this paper we propose a method to compute the solution to the fractional diffusion equation on directed networks, which can be expressed in terms of the graph Laplacian L as a product f (L T) b , where f is a non-analytic function involving fractional powers and b is a given vector. The graph Laplacian is a singular matrix, causing Krylov methods for f (L T) b to converge more slowly. In order to overcome this difficulty and achieve faster convergence, we use rational Krylov methods applied to a desingularized version of the graph Laplacian, obtained with either a rank-one shift or a projection on a subspace. [ABSTRACT FROM AUTHOR]

This paper is concerned with computing approximations of matrix functionals of the form F (A) : = v T f (A) v , where A is a large symmetric positive definite matrix, v is a vector, and f is a Stieltjes function. We approximate F(A) with the aid of rational Gauss quadrature rules. Associated rational Gauss–Radau and rational anti-Gauss rules are developed. Pairs of rational Gauss and rational Gauss–Radau quadrature rules, or pairs of rational Gauss and rational anti-Gauss quadrature rules, can be used to determine upper and lower bounds, or approximate upper and lower bounds, for F(A). The application of rational Gauss rules, instead of standard Gauss rules, is beneficial when the function f has singularities close to the spectrum of A. [ABSTRACT FROM AUTHOR]

Prior recent work, devoted to the study of polynomial Krylov techniques for the approximation of the action of the matrix exponential etAv, is extended to the case of associated φ-functions (which occur within the class of exponential integrators). In particular, a posteriori error bounds and estimates, based on the notion of the defect (residual) of the Krylov approximation are considered. Computable error bounds and estimates are discussed and analyzed. This includes a new error bound which favorably compares to existing error bounds in specific cases. The accuracy of various error bounds is characterized in relation to corresponding Ritz values of A. Ritz values yield properties of the spectrum of A (specific properties are known a priori, e.g., for Hermitian or skew-Hermitian matrices) in relation to the actual starting vector v and can be computed. This gives theoretical results together with criteria to quantify the achieved accuracy on the fly. For other existing error estimates, the reliability and performance are studied by similar techniques. Effects of finite precision (floating point arithmetic) are also taken into account. [ABSTRACT FROM AUTHOR]

Two numerical algorithms are proposed for computing an interval matrix containing the matrix gamma function. In 2014, the author presented algorithms for enclosing all the eigenvalues and basis of invariant subspaces of A ∈ C n × n . As byproducts of these algorithms, we can obtain interval matrices containing diagonal blocks whose spectrums are included in that of A. In this paper, we interpret the interval matrices containing the basis and blocks as a result of verified block diagonalization (VBD), and establish a new framework for enclosing matrix functions using the VBD. To achieve enclosure for the gamma function of the blocks, we derive computable perturbation bounds. We can apply these bounds if disks containing the spectrums of input matrices lie in the open right half plane. We incorporate matrix argument reductions (ARs) to force the input matrices to have this property, and develop theories for accelerating the ARs. The first algorithm uses the VBD based on a numerical spectral decomposition, and involves only cubic complexity if the total computational cost of the accelerated ARs is O (n 3). The second algorithm adopts the VBD based on a numerical Jordan decomposition, is applicable even for defective matrices, and requires O (n 4) operations. [ABSTRACT FROM AUTHOR]

In this paper we are interested in the approximation of fractional powers of self-adjoint positive operators. Starting from the integral representation of the operators, we apply the trapezoidal rule combined with a double-exponential transform of the integrand function. In this work we show how to improve the existing error estimates for the scalar case and also extend the analysis to operators. We report some numerical experiments to show the reliability of the estimates obtained. [ABSTRACT FROM AUTHOR]

The numerical computation of a matrix function such as exp (- t A) V , where A is an n × n large and sparse matrix, V is an n × p block with p ≪ n , and t > 0 arises in various applications including network analysis, the solution of time-dependent partial differential equations (PDE’s) and others. In this work, we propose the use of the global extended-rational Arnoldi method for computing approximations of such functions. The derived method projects the initial problem onto the global extended-rational Krylov subspace RK m e (A , V) = span { ∏ i = 1 m (A + s i I n) - 1 V , … , (A + s 1 I n) - 1 V , V , A V , … , A m - 1 V } of a low dimension. An adaptive procedure of getting the shifts { s 1 , … , s m } during the algorithmic process is given and analyzed. Applications to the solution of time-dependent PDE’s and to network analysis are presented. Numerical examples are presented to show the performance of the global extended-rational Arnoldi process. [ABSTRACT FROM AUTHOR]

Many non-local processes are modeled using mathematical models that include fractional powers of elliptic operators. The approximate solution of stationary problems with fractional powers of operators is most often based on rational approximations introduced in various versions for a fractional power of the self-adjoint positive operator. The purpose of this work is to use such approximations for the approximate solution of nonstationary problems. We consider Cauchy problems for the first and second order differential-operator equations in finite-dimensional Hilbert spaces. Estimates for the proximity of an approximate solution to an exact one are obtained when specifying the absolute and relative errors of the approximation of the fractional power of the operator. We construct splitting schemes based on the additive representation with a rational approximation of the operator's fractional power. The stability and accuracy of factorized two-level additive operator-difference schemes for the first order evolution equation and three-level schemes for a second order equation are established. [ABSTRACT FROM AUTHOR]

In this work, we introduce and investigate a class of matrix-free regularization techniques for discrete linear ill-posed problems based on the approximate computation of a special matrix-function. In order to produce a regularized solution, the proposed strategy employs a regular approximation of the Heavyside step function computed into a small Krylov subspace. This particular feature allows our proposal to be independent from the structure of the underlying matrix. If on the one hand, the use of the Heavyside step function prevents the amplification of the noise by suitably filtering the responsible components of the spectrum of the discretization matrix, on the other hand, it permits the correct reconstruction of the signal inverting the remaining part of the spectrum. Numerical tests on a gallery of standard benchmark problems are included to prove the efficacy of our approach even for problems affected by a high level of noise. [ABSTRACT FROM AUTHOR]

The Krylov subspace method has been investigated and refined for approximating the behaviors of finite or infinite dimensional linear operators. It has been used for approximating eigenvalues, solutions of linear equations, and operator functions acting on vectors. Recently, for time-series data analysis, much attention is being paid to the Krylov subspace method as a viable method for estimating the multiplications of a vector by an unknown linear operator referred to as a transfer operator. In this paper, we investigate a convergence analysis for Krylov subspace methods for estimating operator-vector multiplications. [ABSTRACT FROM AUTHOR]

In a previous paper by the author, a family of iterations for computing the matrix square root was constructed by exploiting a recursion obeyed by Zolotarev's rational minimax approximants of the function z 1 / 2 . The present paper generalizes this construction by deriving rational minimax iterations for the matrix pth root, where p ≥ 2 is an integer. The analysis of these iterations is considerably different from the case p = 2 , owing to the fact that when p > 2 , rational minimax approximants of the function z 1 / p do not obey a recursion. Nevertheless, we show that several of the salient features of the Zolotarev iterations for the matrix square root, including equioscillatory error, order of convergence, and stability, carry over to the case p > 2 . A key role in the analysis is played by the asymptotic behavior of rational minimax approximants on short intervals. Numerical examples are presented to illustrate the predictions of the theory. [ABSTRACT FROM AUTHOR]

The need to compute matrix functions of the form f(A)v, where A ∈ ℝ N × N is a large symmetric matrix, f is a function such that f(A) is well defined, and v≠ 0 is a vector, arises in many applications. This paper is concerned with the situation when A is so large that the evaluation of f(A) is prohibitively expensive. Then, an approximation of f(A)v often is computed by applying a few, say 1 ≤ n ≪ N, steps of the symmetric Lanczos process to A with initial vector v to determine a symmetric tridiagonal matrix T n ∈ ℝ n × n and a matrix V n ∈ ℝ N × n , whose orthonormal columns span a Krylov subspace. The expression Vnf(Tn)e1∥v∥ furnishes an approximation of f(A)v. The evaluation of f(Tn) is inexpensive, because the matrix Tn is small. It is important to be able to estimate the error in the computed approximation. This paper describes a novel approach that is based on a technique proposed by Spalević for estimating the error in Gauss quadrature rules. [ABSTRACT FROM AUTHOR]

Linear and nonlinear evolution equations have been formulated to address problems in various fields of science and technology. Recently, methods using an exponential integrator for solving evolution equations, where matrix functions must be computed repeatedly, have been investigated and refined. In this paper, we propose a new method for computing these matrix functions which is called an inexact rational Krylov method. This is a more efficient version of the rational Krylov method with appropriate shifts, which was proposed by Hashimoto and Nodera (ANZIAM J 58:C149–C161, 2016). The advantage of the inexact rational Krylov method is that it computes linear equations that appear in the rational Krylov method efficiently while guaranteeing the accuracy of the final solution. [ABSTRACT FROM AUTHOR]

In this paper, we study the relative error in the numerical solution of a linear ordinary differential equation y'(t) = Ay(t), t ≥ 0, where A is a normal matrix. The numerical solution is obtained by using at any step an approximation of the matrix exponential, e.g., a polynomial or a rational approximation. The error of the numerical solution with respect to the exact solution is due to this approximation as well as to a possible perturbation in the initial value. For an unperturbed initial value, we have found: (1) unlike the absolute error, the relative error always grows linearly in time; (2) in the long-time, the contributions to the relative error relevant to non-rightmost eigenvalues of A disappear. [ABSTRACT FROM AUTHOR]

Rational Krylov subspaces have become a fundamental ingredient in numerical linear algebra methods associated with reduction strategies. Nonetheless, many structural properties of the reduced matrices in these subspaces are not fully understood. We advance in this analysis by deriving bounds on the entries of rational Krylov reduced matrices and of their functions, that ensure an a-priori decay of their entries as we move away from the main diagonal. As opposed to other decay pattern results in the literature, these properties hold in spite of the lack of any banded structure in the considered matrices. Numerical experiments illustrate the quality of our results. [ABSTRACT FROM AUTHOR]

Evaluating the action of a matrix function on a vector, that is x = f (M) v , is an ubiquitous task in applications. When M is large, one usually relies on Krylov projection methods. In this paper, we provide effective choices for the poles of the rational Krylov method for approximating x when f(z) is either Cauchy–Stieltjes or Laplace–Stieltjes (or, which is equivalent, completely monotonic) and M is a positive definite matrix. Relying on the same tools used to analyze the generic situation, we then focus on the case M = I ⊗ A - B T ⊗ I , and v obtained vectorizing a low-rank matrix; this finds application, for instance, in solving fractional diffusion equation on two-dimensional tensor grids. We see how to leverage tensorized Krylov subspaces to exploit the Kronecker structure and we introduce an error analysis for the numerical approximation of x. Pole selection strategies with explicit convergence bounds are given also in this case. [ABSTRACT FROM AUTHOR]

The backward error analysis is a great tool which allows selecting in an effective way the scaling parameter s and the polynomial degree of approximation m when the action of the matrix exponential exp (A) v has to be approximated by p m (s - 1 A) s v = exp (A + Δ A) v . We propose here a rigorous bound for the relative backward error Δ A 2 / A 2 , which is of particular interest for matrices whose field of values is skinny, such as the discretization of the advection–diffusion or the Schrödinger operators. The numerical results confirm the superiority of the new approach with respect to methods based on the classical power series expansion of the backward error for the matrices of our interest, both in terms of computational cost and achieved accuracy. [ABSTRACT FROM AUTHOR]

Extended Krylov subspace methods are attractive methods for computing approximations of matrix functions and other problems producing large-scale matrices. In this work, we propose the extended nonsymmetric global Lanczos method for solving some matrix approximation problems. The derived algorithm uses short recursive relations to generate bi-orthonormal bases, with respect to the Frobenius inner product, of the corresponding extended Krylov subspaces K m e (A , V) and K m e (A T , W) . Here, A is a large nonsymmetric matrix; V and W ∈ ℝ n × s are two blocks. New algebraic properties of the proposed method are developed and applications to approximation of both WTf(A)V and trace(WTf(A)V) are given. Numerical examples are presented to show the performance of the extended nonsymmetric global Lanczos for these problems. [ABSTRACT FROM AUTHOR]

In this paper, our formulation generalizes the fractional power series to the matrix form and a new version of the matrix fractional Taylor's series is also considered in terms of Caputo's fractional derivative. Moreover, several significant results have been realignment to these generalizations. Finally, to demonstrate the capability and efficiency of our theoretical results, we present the solutions of three linear non-homogenous higher order (m − 1 < α ≤ m, m ∈ N) matrix fractional differential equations by using our new approach. [ABSTRACT FROM AUTHOR]

We study a reliable pole selection for the rational approximation of the resolvent of fractional powers of operators in both the finite and infinite dimensional setting. The analysis exploits the representation in terms of hypergeometric functions of the error of the Padé approximation of the fractional power. We provide quantitatively accurate error estimates that can be used fruitfully for practical computations. We present some numerical examples to corroborate the theoretical results. The behavior of rational Krylov methods based on this theory is also presented. [ABSTRACT FROM AUTHOR]

An a posteriori estimate for the error of a standard Krylov approximation to the matrix exponential is derived. The estimate is based on the defect (residual) of the Krylov approximation and is proven to constitute a rigorous upper bound on the error, in contrast to existing asymptotical approximations. It can be computed economically in the underlying Krylov space. In view of time-stepping applications, assuming that the given matrix is scaled by a time step, it is shown that the bound is asymptotically correct (with an order related to the dimension of the Krylov space) for the time step tending to zero. This means that the deviation of the error estimate from the true error tends to zero faster than the error itself. Furthermore, this result is extended to Krylov approximations of φ -functions and to improved versions of such approximations. The accuracy of the derived bounds is demonstrated by examples and compared with different variants known from the literature, which are also investigated more closely. Alternative error bounds are tested on examples, in particular a version based on the concept of effective order. For the case where the matrix exponential is used in time integration algorithms, a step size selection strategy is proposed and illustrated by experiments. [ABSTRACT FROM AUTHOR]

In this paper, we propose nonoverlapping localized exponential time differencing (ETD) methods for diffusion problems. The model time-dependent diffusion equation is first reformulated on subdomains based on the nonoverlapping domain decomposition, in which Neumann boundary conditions are imposed on the interfaces for the subdomain problems and Dirichlet type conditions are enforced to form a space-time interface problem. After spatial discretization by standard central finite differences and temporal integration with the first or second order ETD methods, the fully discrete interface problem is obtained. Such an interface problem is then solved iteratively either at each time step or over the whole time interval: the former involves the solution of stationary problems in each subdomain at each iteration while the latter involves the solution of time-dependent subdomain problems at each iteration. For both approaches, we prove that localized ETD solutions conserve mass exactly and converge in time to the exact space semidiscrete solution. Numerical experiments in two dimensions are also presented to illustrate the performance of the proposed methods. [ABSTRACT FROM AUTHOR]