Showing papers in "Numerical Algorithms in 2010"
TL;DR: A reduced composition technique has been used on Newton and Jarratt's methods to obtain an optimal relation between convergence order, functional evaluations and number of operations, whose efficiency indices are proved to be better for systems of nonlinear equations.
Abstract: A reduced composition technique has been used on Newton and Jarratt’s methods in order to obtain an optimal relation between convergence order, functional evaluations and number of operations. Following this aim, a family of methods is obtained whose efficiency indices are proved to be better for systems of nonlinear equations.
205 citations
TL;DR: New algorithms and related software for total variation (TV) image reconstruction, more specifically: denoising, inpainting, and deblurring are described, based on one of Nesterov’s first-order methods.
Abstract: This paper describes new algorithms and related software for total variation (TV) image reconstruction, more specifically: denoising, inpainting, and deblurring. The algorithms are based on one of Nesterov’s first-order methods, tailored to the image processing applications in such a way that, except for the mandatory regularization parameter, the user needs not specify any parameters in the algorithms. The software is written in C with interface to Matlab (version 7.5 or later), and we demonstrate its performance and use with examples.
175 citations
TL;DR: An approach to the order conditions based on Butcher’s algebraic theory of trees is presented, and methods that take advantage of cheap computations of the second derivatives are derived.
Abstract: The theory of Runge-Kutta methods for problems of the form y′ = f(y) is extended to include the second derivative y′′ = g(y): = f′(y)f(y) We present an approach to the order conditions based on Butcher’s algebraic theory of trees (Butcher, Math Comp 26:79–106, 1972), and derive methods that take advantage of cheap computations of the second derivatives Only explicit methods are considered here where attention is given to the construction of methods that involve one evaluation of f and many evaluations of g per step Methods with stages up to five and of order up to seven including some embedded pairs are presented The first part of the paper discusses a theoretical formulation used for the derivation of these methods which are also of wider applicability The second part presents experimental results for non-stiff and mildly stiff problems The methods include those with the computation of one second derivative (plus many first derivatives) per step, and embedded methods for changing stepsize as well as those involving one first derivative (plus many second derivatives) per step The experiments have been performed on standard problems and comparisons made with some standard explicit Runge-Kutta methods
125 citations
TL;DR: A family of eighth order methods for the solution of nonlinear equations based on Ostrowski’s fourth order method, which agrees with Kung–Traub conjecture for the case n = 4 and is competitive with other similar robust methods and very effective in high precision computations.
Abstract: Based on Ostrowski’s fourth order method, we derive a family of eighth order methods for the solution of nonlinear equations. In terms of computational cost the family requires three evaluations of the function and one evaluation of first derivative. Therefore, the efficiency index of the present methods is 1.682 which is better than the efficiency index 1.587 of Ostrowski’s method. Kung and Traub conjectured that multipoint iteration methods without memory based on n evaluations have optimal order 2n − 1. Thus, the family agrees with Kung–Traub conjecture for the case n = 4. The efficacy of the present methods is tested on a number of numerical examples. It is observed that our methods are competitive with other similar robust methods and very effective in high precision computations.
103 citations
TL;DR: It can be seen in this paper that the auxiliary parameter $\hbar,$ which controls the convergence of the HAM approximate series solutions, has another important application, predicting and calculating multiple solutions.
Abstract: In this paper, Homotopy Analysis Method (HAM) is applied to numerically approximate the eigenvalues of the fractional Sturm-Liouville problems. The eigenvalues are not unique. These multiple solutions, i.e., eigenvalues, can be calculated by starting the HAM algorithm with one and the same initial guess and linear operator \(\mathcal{L}\). It can be seen in this paper that the auxiliary parameter \(\hbar,\) which controls the convergence of the HAM approximate series solutions, has another important application. This important application is predicting and calculating multiple solutions.
89 citations
TL;DR: Two numerical methods for solving the two-dimensional anomalous subdiffusion equation (2D-ASDE) are presented and their stability, convergence and solvability are discussed.
Abstract: Anomalous dynamics in complex systems have gained much interest in recent years. In this paper, a two-dimensional anomalous subdiffusion equation (2D-ASDE) is considered. Two numerical methods for solving the 2D-ASDE are presented. Their stability, convergence and solvability are discussed. A new multivariate extrapolation is introduced to improve the accuracy. Finally, numerical examples are given to demonstrate the effectiveness of the schemes and confirm the theoretical analysis.
87 citations
TL;DR: It is shown that the orthogonal expansion of the finite-dimensional Wiener process, that appears in the discretized problem, can be truncated severely without losing the asymptotic order of the method, provided that the kernel of the covariance operator of theWiener process is smooth enough.
Abstract: We consider a semilinear parabolic PDE driven by additive noise. The equation is discretized in space by a standard piecewise linear finite element method. We show that the orthogonal expansion of the finite-dimensional Wiener process, that appears in the discretized problem, can be truncated severely without losing the asymptotic order of the method, provided that the kernel of the covariance operator of the Wiener process is smooth enough. For example, if the covariance operator is given by the Gauss kernel, then the number of terms to be kept is the quasi-logarithm of the number of terms in the original expansion. Then one can reduce the size of the corresponding linear algebra problem enormously and hence reduce the computational complexity, which is a key issue when stochastic problems are simulated.
72 citations
TL;DR: The Fréchet derivative of a matrix function f at A in the direction E, where A and E are real matrices, can be approximated by Im f(A + ihE)/h for some suitably small h, and is proved to be of second order in h for analytic functions f and also for the matrix sign function.
Abstract: We show that the Frechet derivative of a matrix function f at A in the direction E, where A and E are real matrices, can be approximated by Im f(A + ihE)/h for some suitably small h. This approximation, requiring a single function evaluation at a complex argument, generalizes the complex step approximation known in the scalar case. The approximation is proved to be of second order in h for analytic functions f and also for the matrix sign function. It is shown that it does not suffer the inherent cancellation that limits the accuracy of finite difference approximations in floating point arithmetic. However, cancellation does nevertheless vitiate the approximation when the underlying method for evaluating f employs complex arithmetic. The ease of implementation of the approximation, and its superiority over finite differences, make it attractive when specialized methods for evaluating the Frechet derivative are not available, and in particular for condition number estimation when used in conjunction with a block 1-norm estimation algorithm.
68 citations
TL;DR: Two criteria for selection are proposed with a view to deriving pairs of all orders 6(5) to 9(8) which minimize computation while achieving a user-specified accuracy.
Abstract: Explicit Runge–Kutta pairs are known to provide efficient solutions to initial value differential equations with inexpensive derivative evaluations. Two criteria for selection are proposed with a view to deriving pairs of all orders 6(5) to 9(8) which minimize computation while achieving a user-specified accuracy. Coefficients of improved pairs, their stability regions and coefficients of appended optimal interpolatory Runge–Kutta formulas are provided on the author’s website (
www.math.sfu.ca/~jverner
). This note reports results of tests on these pairs to illustrate their effectiveness in solving nonstiff initial value problems. These pairs and interpolants may be used for implementation, or else to provide comparison targets for other new types of methods such as explicit general linear methods.
63 citations
TL;DR: This paper employs the Homotopy Analysis Method (HAM) to derive a family of iterative methods for solving systems of nonlinear algebraic equations which are more efficient than their classical counterparts such as Newton's, Chebychev's and Halley's methods.
Abstract: Solving systems of nonlinear equations is a relatively complicated problem for which a number of different approaches have been proposed. In this paper, we employ the Homotopy Analysis Method (HAM) to derive a family of iterative methods for solving systems of nonlinear algebraic equations. Our approach yields second and third order iterative methods which are more efficient than their classical counterparts such as Newton’s, Chebychev’s and Halley’s methods.
59 citations
TL;DR: A new class of two-step Runge-Kutta methods for the numerical solution of ordinary differential equations is proposed by relaxing some of the collocation conditions to obtain methods with desirable stability properties.
Abstract: A new class of two-step Runge-Kutta methods for the numerical solution of ordinary differential equations is proposed. These methods are obtained using the collocation approach by relaxing some of the collocation conditions to obtain methods with desirable stability properties. Local error estimation for these methods is also discussed.
TL;DR: The stability and convergence of an implicit box scheme for subsonic flows modelled by scalar conservation laws with dissipative and possibly stiff source terms is investigated and Kružkov’s theorem is verified.
Abstract: We investigate the stability and convergence of an implicit box scheme for subsonic flows modelled by scalar conservation laws with dissipative and possibly stiff source terms. The scheme is proposed for solving transient gas flow problems in pipeline networks. Such networks are operated in the subsonic flow region and are characterized by pressure losses due to dissipative friction terms. We verify the properties stated by Kružkov’s theorem (Kružkov, Math. USSR-Sb. 10:217–243, 1970) for the approximate solution and prove its convergence to the entropy solution.
TL;DR: An algorithm for computing verified and narrow error bounds with the property that a slightly perturbed system is proved to have a double root within the computed bounds is described.
Abstract: It is well known that it is an ill-posed problem to decide whether a function has a multiple root. Even for a univariate polynomial an arbitrary small perturbation of a polynomial coefficient may change the answer from yes to no. Let a system of nonlinear equations be given. In this paper we describe an algorithm for computing verified and narrow error bounds with the property that a slightly perturbed system is proved to have a double root within the computed bounds. For a univariate nonlinear function f we give a similar method also for a multiple root. A narrow error bound for the perturbation is computed as well. Computational results for systems with up to 1000 unknowns demonstrate the performance of the methods.
TL;DR: Numerical comparisons with conjugate gradient algorithms show that this hybrid computational scheme outperforms a variant of the hybrid conjugates gradient algorithm given by Andrei.
Abstract: An accelerated hybrid conjugate gradient algorithm represents the subject of this paper. The parameter βk is computed as a convex combination of \(\beta_k^{HS}\) (Hestenes and Stiefel, J Res Nat Bur Stand 49:409–436, 1952) and \(\beta_k^{DY}\) (Dai and Yuan, SIAM J Optim 10:177–182, 1999), i.e. \(\beta_k^C =\left({1-\theta_k}\right)\beta_k^{HS} + \theta_k \beta_k^{DY}\). The parameter θk in the convex combinaztion is computed in such a way the direction corresponding to the conjugate gradient algorithm is the best direction we know, i.e. the Newton direction, while the pair (sk, yk) satisfies the modified secant condition given by Li et al. (J Comput Appl Math 202:523–539, 2007) Bk + 1sk = zk, where \(z_k =y_k +\left({{\eta_k} / {\left\| {s_k} \right\|^2}} \right)s_k\), \(\eta_k =2\left( {f_k -f_{k+1}} \right)+\left( {g_k +g_{k+1}} \right)^Ts_k\), sk = xk + 1 − xk and yk = gk + 1 − gk. It is shown that both for uniformly convex functions and for general nonlinear functions the algorithm with strong Wolfe line search is globally convergent. The algorithm uses an acceleration scheme modifying the steplength αk for improving the reduction of the function values along the iterations. Numerical comparisons with conjugate gradient algorithms show that this hybrid computational scheme outperforms a variant of the hybrid conjugate gradient algorithm given by Andrei (Numer Algorithms 47:143–156, 2008), in which the pair (sk, yk) satisfies the classical secant condition Bk + 1sk = yk, as well as some other conjugate gradient algorithms including Hestenes-Stiefel, Dai-Yuan, Polack-Ribiere-Polyak, Liu-Storey, hybrid Dai-Yuan, Gilbert-Nocedal etc. A set of 75 unconstrained optimization problems with 10 different dimensions is being used (Andrei, Adv Model Optim 10:147–161, 2008).
TL;DR: The estimates of the central moments for q-Bernstein operators (0 < q < 1) are given which can be used for studying the approximation properties of the operators.
Abstract: In this note we give the estimates of the central moments for q-Bernstein operators (0 < q < 1) which can be used for studying the approximation properties of the operators.
TL;DR: An algorithm for unconstrained optimization based on the transformation of the Newton method with the line search into a gradient descent method that is proved for uniformly convex functions and strictly convex quadratic functions satisfying specified conditions.
Abstract: We introduced an algorithm for unconstrained optimization based on the transformation of the Newton method with the line search into a gradient descent method. Main idea used in the algorithm construction is approximation of the Hessian by an appropriate diagonal matrix. The steplength calculation algorithm is based on the Taylor’s development in two successive iterative points and the backtracking line search procedure. The linear convergence of the algorithm is proved for uniformly convex functions and strictly convex quadratic functions satisfying specified conditions.
TL;DR: New classes of continuous two-step Runge-Kutta methods for the numerical solution of ordinary differential equations are derived, imposing some interpolation and collocation conditions, in order to obtain desirable stability properties such as A-stability and L-st stability.
Abstract: New classes of continuous two-step Runge-Kutta methods for the numerical solution of ordinary differential equations are derived. These methods are developed imposing some interpolation and collocation conditions, in order to obtain desirable stability properties such as A-stability and L-stability. Particular structures of the stability polynomial are also investigated.
TL;DR: Based on two-grid discretizations, some local and parallel finite element algorithms for the d-dimensional (d = 2,3) transient Stokes equations are proposed and analyzed and errors of the corresponding solutions from these algorithms are estimated.
Abstract: Based on two-grid discretizations, some local and parallel finite element algorithms for the d-dimensional (d = 2,3) transient Stokes equations are proposed and analyzed. Both semi- and fully discrete schemes are considered. With backward Euler scheme for the temporal discretization, the basic idea of the fully discrete finite element algorithms is to approximate the generalized Stokes equations using a coarse grid on the entire domain, then correct the resulted residue using a finer grid on overlapped subdomains by some local and parallel procedures at each time step. By the technical tool of local a priori estimate for the fully discrete finite element solution, errors of the corresponding solutions from these algorithms are estimated. Some numerical results are also given which show that the algorithms are highly efficient.
TL;DR: This paper is the first approach to the solution of Volterra integral equation by exponential fitting methods and develops a Direct Quadrature method, which uses a class of ef-based quadrature rules adapted to the current problem to solve.
Abstract: This paper is the first approach to the solution of Volterra integral equation by exponential fitting methods. We have developed a Direct Quadrature method, which uses a class of ef-based quadrature rules adapted to the current problem to solve. We have analyzed the convergence of the method and have found different formulas for the coefficients, which limit rounding errors for small stepsizes. Numerical experiments for comparison with other DQ methods are presented.
TL;DR: Some recent nodal distributions for a high order interpolation in the triangle are summarized and adopted as approximation points for the numerical solution of an elliptic partial differential equation on an unstructured simplicial mesh.
Abstract: In the field of spectral element approximations, the interpolation points can be chosen on the basis of different criteria, going from the minimization of the Lebesgue constant to the simplicity of the point generation procedure. In the present paper, we summarize some recent nodal distributions for a high order interpolation in the triangle. We then adopt these points as approximation points for the numerical solution of an elliptic partial differential equation on an unstructured simplicial mesh. The L 2-norm of the approximation error is then analyzed for a model problem.
TL;DR: Trees have a central role in the theory of Runge–Kutta methods and they also have applications to more general methods, involving multiple values and multiple stages.
Abstract: This paper presents a review of the role played by trees in the theory of Runge–Kutta methods. The use of trees is in contrast to early publications on numerical methods, in which a deceptively simpler approach was used. This earlier approach is not only non-rigorous, but also incorrect. It is now known, for example, that methods can have different orders when applied to a single equation and when applied to a system of equations; the earlier approach cannot show this. Trees have a central role in the theory of Runge–Kutta methods and they also have applications to more general methods, involving multiple values and multiple stages.
TL;DR: Of notable interest and significance are the uniformity and the simplicity in the SDP formulation that solves effectively many otherwise very difficult QIEPs.
Abstract: In the past decade or so, semi-definite programming (SDP) has emerged as a powerful tool capable of handling a remarkably wide range of problems. This article describes an innovative application of SDP techniques to quadratic inverse eigenvalue problems (QIEPs). The notion of QIEPs is of fundamental importance because its ultimate goal of constructing or updating a vibration system from some observed or desirable dynamical behaviors while respecting some inherent feasibility constraints well suits many engineering applications. Thus far, however, QIEPs have remained challenging both theoretically and computationally due to the great variations of structural constraints that must be addressed. Of notable interest and significance are the uniformity and the simplicity in the SDP formulation that solves effectively many otherwise very difficult QIEPs.
TL;DR: In this paper, the authors proposed two modifications of QN methods based on Newton's and Shamanski's method for singular problems, which belong to the class of two-step iterative rule for matrix updates and the choice of parameters that keep iterative sequence within convergence region.
Abstract: Quasi-Newton methods for solving singular systems of nonlinear equations are considered in this paper. Singular roots cause a number of problems in implementation of iterative methods and in general deteriorate the rate of convergence. We propose two modifications of QN methods based on Newton's and Shamanski's method for singular problems. The proposed algorithms belong to the class of two-step iterations. Influence of iterative rule for matrix updates and the choice of parameters that keep iterative sequence within convergence region are empirically analyzed and some conclusions are obtained.
TL;DR: Under suitable conditions, the modified YT and YT+ methods are proved to prove the global convergence of these two methods.
Abstract: In this paper, by the use of Gram-Schmidt orthogonalization, we propose a class of modified conjugate gradient methods. The methods are modifications of the well-known conjugate gradient methods including the PRP, the HS, the FR and the DY methods. A common property of the modified methods is that the direction generated by any member of the class satisfies \(g_{k}^{T}d_k=-\|g_k\|^2\). Moreover, if line search is exact, the modified method reduces to the standard conjugate gradient method accordingly. In particular, we study the modified YT and YT+ methods. Under suitable conditions, we prove the global convergence of these two methods. Extensive numerical experiments show that the proposed methods are efficient for the test problems from the CUTE library.
TL;DR: In this paper, the Riesz and Moisil-Teodorescu systems were reconstructed in quaternionic language and complete sets of polynomials solutions of both systems that are orthonormal with respect to a certain inner product were constructed.
Abstract: As it is well-known, the generalization of the classical Cauchy-Riemann system to higher dimensions leads to the so-called Riesz and Moisil-Teodorescu systems. Rewriting these systems in quaternionic language and taking advantage of the underlying algebra, we construct complete sets of polynomials solutions of both systems that are orthonormal with respect to a certain inner product. The restrictions of those polynomials to the unit sphere can be viewed as analogues to the complex case of the Fourier exponential functions $\{e^{i n \theta}\}_{n\geq 0}$ on the unit circle and constitute a refinement of the well-known spherical harmonics.
TL;DR: In this paper, the authors developed new methods for solving general linear parabolic partial differential equations (PDEs) in one space dimension, which combine quadratic-spline collocation for space discretization and classical finite differences, such as Crank-Nicolson, for time discretisation.
Abstract: New methods for solving general linear parabolic partial differential equations (PDEs) in one space dimension are developed. The methods combine quadratic-spline collocation for the space discretization and classical finite differences, such as Crank-Nicolson, for the time discretization. The main computational requirements of the most efficient method are the solution of one tridiagonal linear system at each time step, while the resulting errors at the gridpoints and midpoints of the space partition are fourth order. The stability and convergence properties of some of the new methods are analyzed for a model problem. Numerical results demonstrate the stability and accuracy of the methods. Adaptive mesh techniques are introduced in the space dimension, and the resulting method is applied to the American put option pricing problem, giving very competitive results.
TL;DR: Stability domains of these hybrid methods, which are based on backward differentiation formula (BDF) for the numerical solutions of ordinary differential equations (ODEs), are obtained showing that all these new methods, say HBDF, are A(α)-stable whereas they have wide stability regions comparing with those of some known methods.
Abstract: In this article, we have presented the details of hybrid methods which are based on backward differentiation formula (BDF) for the numerical solutions of ordinary differential equations (ODEs). In these hybrid BDF, one additional stage point (or off-step point) has been used in the first derivative of the solution to improve the absolute stability regions. Stability domains of our presented methods have been obtained showing that all these new methods, we say HBDF, of order p, p = 2,4,..., 12, are A(α)-stable whereas they have wide stability regions comparing with those of some known methods like BDF, extended BDF (EBDF), modified EBDF (MEBDF), adaptive EBDF (A-EBDF), and second derivtive Enright methods. Numerical results are also given for five test problems.
TL;DR: A computational algorithm based on the finite element method is introduced and analysed for the solution of a linear first-order functional differential equation which is applicable both in the case of constant and variable coefficients.
Abstract: This paper is devoted to the approximate solution of a linear first-order functional differential equation which involves delayed and advanced arguments. We seek a solution x, defined for t???(0, k???1],(k???IN ), which takes given values on the intervals [???1, 0] and (k???1, k]. Continuing the work started in previous articles on this subject, we introduce and analyse a computational algorithm based on the finite element method for the solution of this problem which is applicable both in the case of constant and variable coefficients. Numerical results are presented and compared with the results obtained by other methods.
TL;DR: The focus here is on solving the minimization subproblem which is part of a method to numerically solve a Monge–Ampère type equation and these subproblems must be evaluated many times in this numerical solution technique and thus efficiency is of utmost importance.
Abstract: This note develops theory and a solution technique for a quadratically constrained eigenvalue minimization problem. This class of problems arises in the numerical solution of fully-nonlinear boundary value problems of Monge–Ampere type. Though it is most important in the three dimensional case, the solution method is directly applicable to systems of arbitrary dimension. The focus here is on solving the minimization subproblem which is part of a method to numerically solve a Monge–Ampere type equation. These subproblems must be evaluated many times in this numerical solution technique and thus efficiency is of utmost importance. A novelty of this minimization algorithm is that it is finite, of complexity \(\mathcal{O}(n^3)\), with the exception of solving a very simple rational function of one variable. This function is essentially the same for any dimension. This result is quite surprising given the nature of the constrained minimization problem.
TL;DR: It is proved that the scheme is second-order convergent with respect to the spatial variable and the Singularity-Separating method is used to remove the singularity of the non-smooth payoff function.
Abstract: In this paper we present a stable numerical method for the linear complementary problem arising from American put option pricing. The numerical method is based on a hybrid finite difference spatial discretization on a piecewise uniform mesh and an implicit time stepping technique. The scheme is stable for arbitrary volatility and arbitrary interest rate. We apply some tricks to derive the error estimates for the direct application of finite difference method to the linear complementary problem. We use the Singularity-Separating method to remove the singularity of the non-smooth payoff function. It is proved that the scheme is second-order convergent with respect to the spatial variable. Numerical results support the theoretical results.