Exact slow-fast decomposition of the singularly perturbed matrix differential Riccati equation

In this paper the Hamiltonian matrix formulation of the Riccati equation is used to derive the reduced-order pure-slow and pure-fast matrix differential Riccati equations of singularly perturbed systems. These pure-slow and pure-fast matrix differential Riccati equations are obtained by decoupling the singularly perturbed matrix differential Riccati equation of dimension n₁+n₂ into the pure-slow regular matrix differential Riccati equation of dimension n₁ and the pure-fast stiff matrix differential Riccati equation of dimension n₂. A formula is derived that produces the solution of the original singularly perturbed matrix differential Riccati equation in terms of solutions of the pure-slow and pure-fast reduced-order matrix differential Riccati equations and solutions of two reduced-order initial value problems. In addition to its theoretical importance, the main result of this paper can also be used to implement optimal filtering and control schemes for singularly perturbed linear time-invariant systems independently in pure-slow and pure-fast time scales.     

A fast Fourier-collocation method for second boundary integral equations

In this paper we develop a fast collocation method for second boundary integral equations by the trigonometric polynomials. We propose a convenient way to compress the dense matrix representation of a compact integral operator with a smooth kernel under the Fourier basis and the corresponding collocation functionals. The compression leads to a sparse matrix with only O(nlog²n) number of nonzero entries, where 2n+1 denotes the order of the matrix. Thus we develop a fast Fourier-collocation method. We prove that the fast Fourier-collocation method gives the optimal convergence order up to a logarithmic factor. Moreover, we design a fast scheme for solving the corresponding truncated linear system. We establish that this algorithm preserves the quasi-optimal convergence of the approximate solution with requiring a number of O(nlog³n) multiplications.     

Fast Construction of Optimal Circulant Preconditioners for Matrices from the Fast Dense Matrix Method

In this paper, we consider solving non-convolution type integral equations by the preconditioned conjugate gradient method. The fast dense matrix method is a fast multiplication scheme that provides a dense discretization matrix A approximating a given integral equation. The dense matrix A can be constructed in O(n) operations and requires only O(n) storage where n is the size of the matrix. Moreover, the matrix-vector multiplication A xcan be done in O(n log n) operations. Thus if the conjugate gradient method is used to solve the discretized system, the cost per iteration is O(n log n) operations. However, for some integral equations, such as the Fredholm integral equations of the first kind, the system will be ill-conditioned and therefore the convergence rate of the method will be slow. In these cases, preconditioning is required to speed up the convergence rate of the method. A good choice of preconditioner is the optimal circulant preconditioner which is the minimizer of C – A _F in Frobenius norm over all circulant matrices C. It can be obtained by taking arithmetic averages of all the entries of A and therefore the cost of constructing the preconditioner is of O(n ²) operations for general dense matrices. In this paper, we develop an O(n log n) method of constructing the preconditioner for dense matrices A obtained from the fast dense matrix method. Application of these ideas to boundary integral equations from potential theory will be given. These equations are ill-conditioned whereas their optimal circulant preconditioned equations will be well-conditioned. The accuracy of the approximation A, the fast construction of the preconditioner and the fast convergence of the preconditioned systems will be illustrated by numerical examples.     

Matrix elements of vertex operators of the deformed W A n -algebra and the Harish-Chandra solutions to Macdonald's difference equations

In this paper we prove that certain matrix elements of vertex operators of the deformed W A _n -algebra satisfy Macdonald's difference equations and form a natural (n + 1)!-dimensional space of solutions. These solutions are the analogues of the Harish-Chandra solutions of the radial parts of the Laplace-Casimir operators on noncompact Riemannian symmetric spaces G/K with prescribed asymptotic behavior. We obtain formulas for analytic continuation of our Harish-Chandra type solutions as a consequence of braiding properties (obtained earlier by Y. Asai, M. Jimbo, T. Miwa, and Y. Pugay) of certain vertex operators of the deformed W A _n -algebra.     

Numerical investigation into the dependence of the Allen–Cahn equation on the free energy

In this paper, we study a direct parallel-in-time (PinT) algorithm for first- and second-order time-dependent differential equations. We use a second-order boundary value method as the time integrator. Instead of solving the corresponding all-at-once system iteratively, we diagonalize the time discretization matrix B, which yields a direct parallel implementation across all time steps. A crucial issue of this methodology is how the condition number (denoted by Cond₂(V )) of the eigenvector matrix V of B behaves as n grows, where n is the number of time steps. A large condition number leads to large roundoff error in the diagonalization procedure, which could seriously pollute the numerical accuracy. Based on a novel connection between the characteristic equation and the Chebyshev polynomials, we present explicit formulas for V and V^− 1, by which we prove that Cond\(_{2}(V)=\mathcal {O}(n^{2})\). This implies that the diagonalization process is well-conditioned and the roundoff error only increases moderately as n grows, and thus, compared to other direct PinT algorithms, a much larger n can be used to yield satisfactory parallelism. A fast structure-exploiting algorithm is also designed for computing the spectral diagonalization of B. Numerical results on parallel machine are given to support our findings, where over 60 times speedup is achieved with 256 cores.

     

An extended fast algorithm for constructing the Dixon resultant matrix

In recent years,the Dixon resultant matrix has been used widely in the re-sultant elimination to solve nonlinear polynomial equations and many researchers havestudied its efficient algorithms.The recursive algorithm is a very efficient algorithm,butwhich deals with the case of three polynomial equations with two variables at most.Inthis paper,we extend the algorithm to the general case of n 1 polynomial equations in nvariables.The algorithm has been implemented in Maple 9.By testing the random polyno-mial equations,the results demonstrate that the efficiency of our program is much betterthan the previous methods,and it is exciting that the necessary condition for the existenceof common intersection points on four general surfaces in which the degree with respectto every variable is not greater than 2 is given out in 48×48 Dixon matrix firstly by ourprogram.     

Parallel Dixon Matrices by Bracket

It is known that the Dixon matrix can be constructed in parallel either by entry or by diagonal. This paper presents another parallel matrix construction, this time by bracket. The parallel by bracket algorithm is the fastest among the three, but not surprisingly it requires the highest number of processors. The method also shows analytically that the Dixon matrix has a total of m(m+1)²(m+2)n(n+1)²(n+2)/36 brackets but only mn(m+1)(n+1)(mn+2m+2n+1)/6 of them are distinct.     

Asynchronous Multisplitting GAOR Method and Asynchronous Multisplitting SSOR Method for Systems of Weakly Nonlinear Equations

In this article, we introduce two new asynchronous multisplitting methods for solving the system of weakly nonlinear equations Ax = G(x) in which A is an n × n real matrix and G(x) = (g ₁(x), g ₂(x), . . . , g _n (x)) ^T is a P-bounded mapping. First, by generalized accelerated overrelaxation (GAOR) technique, we introduce the asynchronous parallel multisplitting GAOR method (including the synchronous parallel multisplitting AOR method as a special case) for solving the system of weakly nonlinear equations. Second, asynchronous parallel multisplitting method based on symmetric successive overrelaxation (SSOR) multisplitting is introduced, which is called asynchronous parallel multisplitting SSOR method. Then under suitable conditions, we establish the convergence of the two introduced methods. The given results contain synchronous multisplitting iterations as a special case.     

Optimal Control and Geodesics on Quadratic Matrix Lie Groups

The purpose of this paper is to extend the symmetric representation of the rigid body equations from the group SO (n) to other groups. These groups are matrix subgroups of the general linear group that are defined by a quadratic matrix identity. Their corresponding Lie algebras include several classical semisimple matrix Lie algebras. The approach is to start with an optimal control problem on these groups that generates geodesics for a left-invariant metric. Earlier work by Bloch, Crouch, Marsden, and Ratiu defines the symmetric representation of the rigid body equations, which is obtained by solving the same optimal control problem in the particular case of the Lie group SO (n). This paper generalizes this symmetric representation to a wider class of matrix groups satisfying a certain quadratic matrix identity. We consider the relationship between this symmetric representation of the generalized rigid body equations and the generalized rigid body equations themselves. A discretization of this symmetric representation is constructed making use of the symmetry, which in turn give rise to numerical algorithms to integrate the generalized rigid body equations for the given class of matrix Lie groups. Dedicated to Professor Arieh Iserles on the Occasion of his Sixtieth Birthday.     

Inner-Outer Factorization via Lur'e Equations

In this paper we construct inner-outer factorizations of transfer functions governed by linear time-invariant differential-algebraic systems. This construction is based on the solution of certain Lur’e equations. In contrast to previous work we do not assume any condition apart from behavioral stabilizability of the underlying system realization. (© 2016 Wiley-VCH Verlag GmbH & Co. KGaA, Weinheim)     

A fast Fourier-Galerkin method for solving integral equations of second kind with weakly singular kernels

We propose in this paper a convenient way to compress the dense matrix representation of a compact integral operator with a weakly singular kernel under the Fourier basis. This compression leads to a sparse matrix with only ${\mathcal{O}}(n\log n)$ number of nonzero entries, where 2n+1 denotes the order of the matrix. Based on this compression strategy, we develop a fast Fourier-Galerkin method for solving second kind integral equations with weakly singular kernels. We prove that the approximate order of the truncated equation remains optimal and that the spectral condition number of the coefficient matrix of the truncated linear system is uniformly bounded. Furthermore, we develop a fast algorithm for solving the corresponding truncated linear system, which preserves the optimal order of the approximate solution with only ${\mathcal{O}}(n\log^{2}n)$ number of multiplications required. Numerical examples complete the paper.     

On computing minimal realizable spectral radii of non‐negative matrices

For decades considerable efforts have been exerted to resolve the inverse eigenvalue problem for non‐negative matrices. Yet fundamental issues such as the theory of existence and the practice of computation remain open. Recently, it has been proved that, given an arbitrary (n–1)‐tuple ?? = (λ₂,…,λ_n) ∈ ?^n–1 whose components are closed under complex conjugation, there exists a unique positive real number ?(??), called the minimal realizable spectral radius of ??, such that the set {λ₁,…,λ_n} is precisely the spectrum of a certain n × n non‐negative matrix with λ₁ as its spectral radius if and only if λ₁ ? ?(??). Employing any existing necessary conditions as a mode of checking criteria, this paper proposes a simple bisection procedure to approximate the location of ?(??). As an immediate application, it offers a quick numerical way to check whether a given n‐tuple could be the spectrum of a certain non‐negative matrix. Copyright © 2004 John Wiley & Sons, Ltd.     

Reduktionsverfahren für Differenzengleichungen bei Randwertaufgaben II

Summary In a recent paper [11], two of the authors investigated a fast reduction method for solving difference equations which approximate certain boundary value problems for Poisson's equation. In this second paper, we prove the numerical stability of the reduction method, and also report on further developments of the method. For the general case, the provided bounds for the numerical errors behave roughly like the condition numberO(n ²) of the linear system; for more realistic model problems estimates of order less thanO(n) are obtained (n ^–1=h=mesh width). The number of operations required for the reduction method isO(n ² ), for the usual five-point difference formula, as well as for the common nine-point formula with discretization error of orderh ⁴.     

On splitting of linear nonstationary systems of differential equations

Using the matrix asymptotic method, we solve the problem of the splitting of n-order linear nonstationary systems of differential equations into n first-order independent equations for the case where the eigenvalues m_i( t) 1 0,i = 1,...,n,t ? [ 0,l ] {\mu_i}\left( \tau \right) \ne 0,i = 1,...,n,\tau \in \left[ {0,\ell } \right] , of the matrix U( t) U\left( \tau \right) determined by the initial system of differential equations are identically equal.     

Finding all solutions to polynomial systems and other systems of equations

In a previous paper, the authors suggested a procedure for obtaining all solutions to certain systems ofn equations inn complex variables. The idea was to start with a trivial system of equations to which all solutions were easily known. The trivial system was then perturbed into the given system. During the perturbation process, one followed the solution paths from each of the trivial solutions into the solutions of the given system. All solutions to the given system were thereby obtained.This paper utilizes a different approach that eliminates the requirement of the previous paper for a leading dominating term in each equation. We add a dominating term artificially and then fade it. Also we rely on mathematically more fundamental concepts from differential topology. These advancements permit the calculation of all solutions to arbitrary polynomials and to various other systems ofn equations inn complex variables. In addition, information on the number of solutions can be obtained without calculation.Work supported in part by NSF Grant No. MCS77-15509 and ARO Grant No. DAAG-29-78-G-0160.Work supported in part by ARO Grant No. DAAG-29-78-G-0160     

On richardson extrapolation for finite difference methods on regular grids

Summary Difference solutions of partial differential equations can in certain cases be expanded by even powers of a discretization parameterh. If we haven solutions corresponding to different mesh widthsh ₁,...,h _n we can improve the accuracy by Richardson extrapolation and get a solution of order 2n, yet only on the intersection of all grids used, i.e. normally on the coarsest grid. To interpolate this high order solution with the same accuracy in points not belonging to all grids, we need 2n points in an interval of length (2n–1)h ₁.This drawback can be avoided by combining such an interpolation with the extrapolation byh. In this case the approximation depends only on grid points in an interval of length 3/2h ₁. The length of this interval is independent of the desired order.By combining this approach with the method of Kreiss, boundary conditions on curved boundaries can be discretized with a high order even on coarse grids.This paper is based on a lecture held at the 5th Sanmarinian University Session of the International Academy of Sciences San Marino, at San Marino, 1988-08-27-1988-09-05     

Intermediate-thrust arcs and their optimality in a central,time-invariant force field

This paper presents the general equations of the intermediate-thrust arcs in a general, time-invariant, central force field. Two families of planar arcs, namely, the family of Lawden's spirals in the equatorial plane of an oblate planet and the family of intermediate-thrust arcs in a gravitational field of the form /r ⁿ , have been considered in detail. The Kelley-Contensou condition has been used to test their optimality condition. It is shown that, in the first case, there exist portions of the arcs at a finite distance satisfying the condition, while, in the second case, the entire family satisfies the condition forn 3. Hence, in a perturbed Newtonian gravitational force field, the intermediate-thrust arcs, under certain favorable conditions, can be part of an optimal trajectory.     

Reduced Vertex Set Result for Interval Semidefinite Optimization Problems

In this paper we propose a reduced vertex result for the robust solution of uncertain semidefinite optimization problems subject to interval uncertainty. If the number of decision variables is m and the size of the coefficient matrices in the linear matrix inequality constraints is n×n, a direct vertex approach would require satisfaction of 2 ^{n(m+1)(n+1)/2} vertex constraints: a huge number, even for small values of n and m. The conditions derived here are instead based on the introduction of m slack variables and a subset of vertex coefficient matrices of cardinality 2 ⁿ⁻¹, thus reducing the problem to a practically manageable size, at least for small n. A similar size reduction is also obtained for a class of problems with affinely dependent interval uncertainty. This work is supported by MIUR under the FIRB project “Learning, Randomization and Guaranteed Predictive Inference for Complex Uncertain Systems,” and by CNR RSTL funds.     

Method of dual matrices for function minimization

In this paper, the method of dual matrices for the minimization of functions is introduced. The method, which is developed on the model of a quadratic function, is characterized by two matrices at each iteration. One matrix is such that a linearly independent set of directions can be generated, regardless of the stepsize employed. The other matrix is such that, at the point where the first matrix fails to yield a gradient linearly independent of all the previous gradients, it generates a displacement leading to the minimal point. Thus, the one-dimensional search is bypassed. For a quadratic function, it is proved that the minimal point is obtained in at mostn + 1 iterations, wheren is the number of variables in the function. Since the one-dimensional search is not needed, the total number of gradient evaluations for convergence is at mostn + 2.Three algorithms of the method are presented. A reverse algorithm, which permits the use of only one matrix, is also given. Considerations pertaining to the applications of this method to the minimization of a quadratic function and a nonquadratic function are given. It is believed that, since the one-dimensional search can be bypassed, a considerable amount of computational saving can be achieved.This paper, supported by the National Science Foundation, Grant No. GP-32453, is based on Ref. 1.     

Halley and Newton are one step apart

One of the central problems of scientific computation is the efficient numerical solution of the system of n equations in n unknowns F (x) = 0 where F: Rⁿ → Rⁿ is sufficiently smooth. While Newton's method is usually used for solving such systems, third order methods will in general use fewer iterations than a second order method to reach the same accuracy. However, the number of arithmetic operations per iteration is higher for third order methods than second order methods. In this note we will consider the case where F = ∇f, where f is three times continuously differentiable. We will show that for a large class of sparse problems the ratio of the number of arithmetic operations of a third order method and Newton's method is constant per iteration. It is shown that when the structure of the tensor is induced by a general sparse structured Hessian matrix which gives no fill-ins when we use a direct method to solve a system of linear equations. (© 2008 WILEY-VCH Verlag GmbH & Co. KGaA, Weinheim)