Discrete Stochastic Arithmetic for Validating Results of Numerical Software

The Discrete Stochastic Arithmetic DSA is a probabilistic approach for round-off error propagation. After a brief review of the CESTAC (Controle et Estimation Stochastique des Arrondis de Calculs) method, which is the basis of DSA, the concept of the “informatical zero”, also called “computational zero”, is defined. The stochastic order relations of the DSA are presented. The DSA is the joint use of the synchronous implementation of the CESTAC method and the stochastic order relations. After having summarized the asynchronous implementation of the CESTAC method, which has been used in the Prosolver software, and which has been legitimately criticized, the synchronous implementation is presented. Then the CADNA (Control of Accuracy and Debugging for Numerical Application) library which implements the DSA arithmetic is presented. It is shown that this library is able to dynamically control the validity of the hypotheses which must hold so that results provided by CESTAC method are reliable. If the hypotheses do not hold then warnings are printed in a special file. The user is informed that numerical anomalies have been detected. Depending on these warnings the user may conclude either that the results obtained are not reliable and that they cannot be correctly computed with this computer, or he may try to debug his code. It is shown that the numerical examples that support the criticisms and which make the Prosolver software fail, do not jeopardize the CADNA library.     

Reducing and monitoring round-off error propagation for symplectic implicit Runge-Kutta schemes

We propose an implementation of symplectic implicit Runge-Kutta schemes for highly accurate numerical integration of non-stiff Hamiltonian systems based on fixed point iteration. Provided that the computations are done in a given floating point arithmetic, the precision of the results is limited by round-off error propagation. We claim that our implementation with fixed point iteration is near-optimal with respect to round-off error propagation under the assumption that the function that evaluates the right-hand side of the differential equations is implemented with machine numbers (of the prescribed floating point arithmetic) as input and output. In addition, we present a simple procedure to estimate the round-off error propagation by means of a slightly less precise second numerical integration. Some numerical experiments are reported to illustrate the round-off error propagation properties of the proposed implementation.     

The equality relations in scientific computing

The equality relation (more generally, the ordering relations) in floating point arithmetic is the exact translation of the mathematical equality relation. Because of the propagation of round-off errors, the floating point arithmetic is not the exact representation of the theoretical arithmetic which is continuous on the real numbers.This leads to some incoherence when the equality concept is used in floating point arithmetic. A well known example is the detection of a zero element in the pivoting column and equation when applying Gaussian elimination, which is almost impossible in floating point arithmetic.We shall begin by showing the inadequacy of the equality relation used in floating point arithmetic (we will call it floating point equality), and then introduce two new concepts: stochastic numbers and the equality relation between such numbers which will be called the stochastic equality. We will show how these concepts allow to recover the coherence between the arithmetic operators and the ordering relations that was missing in floating point computations.     

A posteriori error analysis of round-off errors in the numerical solution of ordinary differential equations

We prove sharp, computable error estimates for the propagation of errors in the numerical solution of ordinary differential equations. The new estimates extend previous estimates of the influence of data errors and discretization errors with a new term accounting for the propagation of numerical round-off errors, showing that the accumulated round-off error is inversely proportional to the square root of the step size. As a consequence, the numeric precision eventually sets the limit for the pointwise computability of accurate solutions of any ODE. The theoretical results are supported by numerically computed solutions and error estimates for the Lorenz system and the van der Pol oscillator.     

Propagation of round-off errors and the role of stability in numerical methods for linear and nonlinear PDEs

It is shown that the customary assumption on the propagation of round-off errors in numerical methods for PDEs is unrealistic, as it yields a convergence result which is better than the best possible similar convergence result for ODEs. A solution is suggested by which round-off errors can be modelled by smooth functions, with consequent weakening of overall stability conditions and improvement of convergence conditions.     

Intrinsic Compiler Support for Interval Arithmetic

Interval arithmetic provides an efficient method for monitoring errors in numerical computations and for solving problems that cannot be efficiently solved with floating-point arithmetic. To support interval arithmetic, several software tools have been developed including interval arithmetic libraries, extended scientific programming languages, and interval-enhanced compilers. The main disadvantage of these software tools is their speed, since interval operations are implemented using function calls. In this paper, compiler support for interval arithmetic is investigated. In particular, the performance benefits of having the compiler inline interval operations to eliminate function call overhead is researched. Interval operations are inlined with the GNU gcc compiler and the performance of interval arithmetic is evaluated on a superscalar architecture. To implement interval operations with compiler support, the compiler produces sequences of instructions that use existing floating point hardware. Simulation results show that the compiler implementation of interval arithmetic is approximately 4 to 5 times faster than a functionally equivalent interval arithmetic software implementation with function call overhead and approximately 1.2 to 1.5 times slower than a dedicated interval arithmetic hardware implementation.     

General algorithm of computation of <Emphasis Type="Italic">c</Emphasis>-table and detection of valleys

We present a review of all interesting results concerning the c-table obtained by the authors for the last two decades. These results are not widely known because they were presented in publications of limited circulation. We discuss different computational aspects of software producing the c-tables in the presence of blocs and their evolution following the evolution of the computer environment: effects of the use of 32-bit arithmetic .≈8 digits), 64-bit arithmetic (double precision, ≈16 digits), and Bailey’s Fortran multiprecision package .32 or 64 digits), competition between the ascending and descending algorithms, relationship between the complexity of computation and precision, overflow and underflow problems, competition between different formulas allowing one to overcome the blocs in the c-table, practical simple criterion of detecting numerical zeros in the c-table allowing to identify the blocs, and automatic detection of valleys.     

回归型支持向量机调节熵函数法的区间扩张研究

为减小由于二进制编码的舍入误差对该问题计算结果的影响,对求解回归支持向量机的一种调节熵方法进行了区间扩张,讨论了区间函数的相关定理与收敛性.对设计的区间算法做了收敛性证明,并给出了数值实验,验证了方法与算法的可行性和有效性.     

Error analysis in stochastic solutions of population balance equations

Stochastic simulation of population balance equations (PBEs) is robust and flexible; however, it exhibits intrinsic stochastic errors which decreases at a very slow rate when increasing the computational resolution. Generally, these stochastic methods can be classified into two groups: (i) the classical Gillespie method and (ii) weighted flow algorithm. An analytical relationship is derived for the first time to connect the variances in these two groups. It also provides a detailed analysis of the resampling process, which has not been given appropriate attention previously. It is found that resampling has a profound effect on the numerical precision. Moreover, by comparing the time evolutions between systematic errors (i.e., errors in the mean value) and stochastic errors (i.e., variances), it is found that the former grows considerably faster than the latter; thus, systematic errors eventually dominate. The present findings facilitate the choice of the most suitable stochastic method for a specific PBE a priori in order to balance numerical precision and efficiency.     

Adaptive Precision Floating-Point Arithmetic and Fast Robust Geometric Predicates

Exact computer arithmetic has a variety of uses, including the robust implementation of geometric algorithms. This article has three purposes. The first is to offer fast software-level algorithms for exact addition and multiplication of arbitrary precision floating-point values. The second is to propose a technique for adaptive precision arithmetic that can often speed these algorithms when they are used to perform multiprecision calculations that do not always require exact arithmetic, but must satisfy some error bound. The third is to use these techniques to develop implementations of several common geometric calculations whose required degree of accuracy depends on their inputs. These robust geometric predicates are adaptive; their running time depends on the degree of uncertainty of the result, and is usually small. These algorithms work on computers whose floating-point arithmetic uses radix two and exact rounding, including machines complying with the IEEE 754 standard. The inputs to the predicates may be arbitrary single or double precision floating-point numbers. C code is publicly available for the two-dimensional and three-dimensional orientation and incircle tests, and robust Delaunay triangulation using these tests. Timings of the implementations demonstrate their effectiveness. Received May 16, 1996, and in revised form March 10, 1997.     

A stochastic approach to global error estimation in ODE multistep numerical integration

In order to assess the quality of approximate solutions obtained in the numerical integration of ordinary differential equations related to initial-value problems, there are available procedures which lead to deterministic estimates of global errors. The aim of this paper is to propose a stochastic approach to estimate the global errors, especially in the situations of integration which are often met in flight mechanics and control problems. Treating the global errors in terms of their orders of magnitude, the proposed procedure models the errors through the distribution of zero-mean random variables belonging to stochastic sequences, which take into account the influence of both local truncation and round-off errors. The dispersions of these random variables, in terms of their variances, are assumed to give an estimation of the errors. The error estimation procedure is developed for Adams-Bashforth-Moulton type of multistep methods. The computational effort in integrating the variational equations to propagate the error covariance matrix associated with error magnitudes and correlations is minimized by employing a low-order (first or second) Euler method. The diagonal variances of the covariance matrix, derived using the stochastic approach developed in this paper, are found to furnish reasonably precise measures of the orders of magnitude of accumulated global errors in short-term as well as long-term orbit propagations.     

A finite element lumped mass scheme for solving eigenvalue problems of circular arches

Summary This study is a continuation of a previous paper [4] in which the numerical results are given by using single precision arithmetic. In this paper, we show the numerical results which experess the sharper convergence properties than those of [4], by using double precision arithmetic.Dedicated to Prof. Masaya Yamaguti on the occasion of his 60th birthday     

Analysis of the GCR method with mixed precision arithmetic using QuPAT

To verify computation results of double precision arithmetic, a high precision arithmetic environment is needed. However, it is difficult to use high precision arithmetic in ordinary computing environments without any special hardware or libraries. Hence, we designed the quadruple precision arithmetic environment QuPAT on Scilab to satisfy the following requirements: (i) to enable programs to be written simply using quadruple precision arithmetic; (ii) to enable the use of both double and quadruple precision arithmetic at the same time; (iii) to be independent of any hardware and operating systems.To confirm the effectiveness of QuPAT, we applied the GCR method for ill-conditioned matrices and focused on the scalar parameters α and β in GCR, partially using DD arithmetic. We found that the use of DD arithmetic only for β leads to almost the same results as when DD arithmetic is used for all computations. We conclude that QuPAT is an excellent interactive tool for using double precision and DD arithmetic at the same time.     

Round-off error for retarded ordinary differential equations:A priori bounds and estimates

Summary Consider the numerical solution of a retarded ordinary differential equation (RODE) by some standard algorithms. For a linear RODE, we estimate the accumulated round-off error as a linear combination of the preceding local round-off errors, and we bound the accumulated round-off error. For a non-linear RODE, we obtain by linearization similar estimates and bounds for the dominant part of the accumulated round-off error.Presented at SIAM National Meeting, June, 1971, Seattle, Washington.     

Software to support numerical analysis teaching

This paper describes a Fortran90 library designed to support the teaching of numerical analysis and its applications. As well as covering traditional material it introduces recent and important ideas in numerical computation such as interval arithmetic and automatic differentiation. The library rests on a module realpac which provides real arithmetic in a range of precisions with a choice of rounding strategies. This, in turn, supports the implementation of an interval arithmetic module intpac. Derived data types and overloaded operations help inexperienced users to interface with unfamiliar data types such as intervals. The library also includes more conventional modules such as lepac for solving linear systems and minpac for nonlinear optimization. These, however, can be enhanced by being linked to more sophisticated tools for sparse matrix handling and automatic differentiation. As well as showing the main structure and scope of the software, the paper mentions some exercises that have successfully been performed by students.     

Stochastic approach of the permutation-perturbation method for round-off error analysis

Any result of algorithms provided by a computer always contains an error resulting from floating-point arithmetic round-off error propagation. The Permutation-Perturbation method is a very efficient practical method for evaluating this error and consequently for estimating the exact significant decimal figures of any result. The stochastic approach of this method, its probabilistic proof, and the perfect agreement between the theoretical and practical aspects are presented in this paper.     

Calculation of Stresses and Consistent Tangent Moduli from Automatic Differentiation of Hyerelastic Strain Energy Functions Through the Use of Hyper Dual Numbers

Many materials as e.g. engineering rubbers, polymers and soft biological tissues are often described by hyperelastic strain energy functions. For their finite element implementation the stresses and consistent tangent moduli are required and obtained mainly in terms of the first and second derivative of the strain energy function. Depending on its mathematical complexity in particular for anisotropic media the analytic derivatives may be expensive to be calculated or implemented. Then numerical approaches may be a useful alternative reducing the development time. Often-used classical finite difference schemes are however quite sensitive with respect to perturbation values and they result in a poor accuracy. The complex-step derivative approximation does never suffer from round-off errors, cf. [1], [2], but it can only provide first derivatives. A method which also provides higher order derivatives is based on hyper dual numbers [3]. This method is independent on the choice of perturbation values and does thus neither suffer from round-off errors nor from approximation errors. Therefore, here we make use of hyper dual numbers and propose a numerical scheme for the calculation of stresses and tangent moduli which are almost identical to the analytic ones. Its uncomplicated implementation and accuracy is illustrated by some representative numerical examples. (© 2014 Wiley-VCH Verlag GmbH & Co. KGaA, Weinheim)     

赌博破产概率及其随机模拟试验

讨论了赌博问题中的最终破产概率,并给出了破产概率的随机模拟计算流程和一个具体例子的数值模拟结果.计算结果表明,由此方法得到破产概率的估计值与理论值的误差很小.最后,通过随机模拟给出游戏结束的平均次数.     

Minimal condition number for positive definite Hankel matrices using semidefinite programming

We present a semidefinite programming approach for computing optimally conditioned positive definite Hankel matrices of order n. Unlike previous approaches, our method is guaranteed to find an optimally conditioned positive definite Hankel matrix within any desired tolerance. Since the condition number of such matrices grows exponentially with n, this is a very good test problem for checking the numerical accuracy of semidefinite programming solvers. Our tests show that semidefinite programming solvers using fixed double precision arithmetic are not able to solve problems with n>30. Moreover, the accuracy of the results for 24?n?30 is questionable. In order to accurately compute minimal condition number positive definite Hankel matrices of higher order, we use a Mathematica 6.0 implementation of the SDPHA solver that performs the numerical calculations in arbitrary precision arithmetic. By using this code, we have validated the results obtained by standard codes for n?24, and we have found optimally conditioned positive definite Hankel matrices up to n=100.     

Estimating manipulation errors in finite element analysis part II

With the increasing use of microcomputers for engineering research and production analyses, control of round-off error is of increasing concern. An essential ingredient in control is the ability to estimate the errors. Part I of this paper described the errors in processing numbers with digital computers [4]. This part, Part II, presents and illustrates the effectiveness of error estimation for three goals: pre-estimating, measuring actual errors, and identifying the error source. Pre-estimation is intended to eliminate analyses whose results would be inaccurate. Measurement of actual errors qualifies the accuracy of calculations. Error source identification leads to selection of error reduction methods.The authors conclude that bound estimates of lost precision are excessively high, do not correlate with actual round-off errors, and are not reliable in identifying the key sources of error. Experiments show that the closure measure provides excellent estimates of maximum error in stress predictions. This fact suggests a conservative approach to protect against use of inaccurate analysis results.