Minimax Polynomial Approximation Research Articles

Error analyses and evaluation for solving burnup equations with MMPA method

The Mini-Max Polynomial Approximation (MMPA) is a new method to solve burnup equations, which has attracted much attention due to its advantages of good computational accuracy and efficiency. Understanding its error mechanism is significant for its application in practical burnup calculations. In this research, a mathematical model is proposed to clarify the error mechanism of the MMPA method for the burnup calculation. Then the coefficients of the MMPA method with different orders are given by the Remez algorithm. The errors of the MMPA method with different orders in burnup calculation are compared and analyzed, especially in single-step burnup calculation and infinite-step burnup calculation. Three cases are selected to verify the correctness of the mathematical error model, including fresh fuel burnup, depleted fuel burnup and decay calculation. The results show that when the order is greater than 36, the calculation accuracy of the MMPA method is less than the numerical noise of double precision calculation. As the order increases, the step length sensitivity of the MMPA method decreases. Finally, the computational accuracy and computational efficiency of the MMPA method with different orders are compared with those of the 16th order Chebyshev Rational Approximation Method (CRAM). The accuracy of the 16th order CRAM is between the 32nd and 36th order MMPA method. And the computation time of MMPA methods below order 76 is less than that of the 16th order CRAM. The theoretical analysis results of the mathematical error model are consistent with the numerical results. Then several best practices for solving burnup equations using the MMPA method are presented: (1) It is recommended to use the MMPA method with orders above 36 and below 76, whose accuracy and efficiency are better than 16th order CRAM in different degrees. (2) In multi-physics coupled computation which may require a strategy of large step length calculation, increasing the order of the MMPA method can not only improve the computational accuracy but also improve the computational efficiency.

Investigation and improvement of the Mini-Max Polynomial Approximation method for solving burnup equations

In this work, a new efficient solution algorithm named Combined Matrix Iterative Algorithm based on the Mini-Max Polynomial Approximation (MMPA) method is proposed for burnup calculations. As a new matrix exponential method in recent years, MMPA method has good numerical stability and computational accuracy, and all calculation are real number operations. As the same with Chebyshev Rational Approximation Method (CRAM) which is widely and successfully used in solving burnup equations, the method can also directly deal with rigid and sparse burnup matrices. Compared with the traditional solution algorithms based on MMPA method, the Combined Matrix Iterative Algorithm (MMPA-CMIA) proposed in this paper can effectively reduce the computational complexity and improve the computational efficiency after equivalently transforming the expression of the MMPA method. Then the MMPA-CMIA is used to theoretically demonstrate and analyze the solution efficiency of random matrices with different orders which has the similar distribution characteristics to burnup matrices. Compared with the 16-order CRAM method and the other two traditional solution algorithms based on MMPA method, the algorithm shows good acceleration effect and computational precision. Finally, the algorithm is used in the actual burnup benchmarks to verify and validate the accuracy and efficiency. The burnup matrices of different orders are constructed through detailed burnup chain and simplified burnup chain. The results show that the algorithm has good accuracy and numerical stability. And comparing with 16-order CRAM, the algorithm can save about 30% of computing time for solving burnup equations. It shows that the algorithm is a potentially powerful method for accelerating the solution of large-scale burnup calculations.

Solving Burnup Equations by Numerical Inversion of the Laplace Transform Using Padé Rational Approximation

Burnup calculations play a very important role in nuclear reactor design and analysis, and solving burnup equations is an essential topic in burnup calculations. In the last decade, several high-accuracy methods, mainly including the Chebyshev rational approximation method (CRAM), the quadrature-based rational approximation method, the Laguerre polynomial approximation method, and the mini-max polynomial approximation method, have been proposed to solve the burnup equations. Although these methods have been demonstrated to be quite successful in the burnup calculations, limitations still exist in some cases, one of which is that the accuracy becomes compromised when treating the time-dependent polynomial external feed rate. In this work, a new method called the Padé rational approximation method (PRAM) is proposed. Without directly approximating the matrix exponential, this new method is derived by using the Padé rational function to approximate the scalar exponential function in the formula of the inverse Laplace transform of burnup equations. Several test cases are carried out to verify the proposed new method. The high accuracy of the PRAM is validated by comparing the numerical results with the high-precision reference solutions. Against CRAM, PRAM is significantly superior in handling the burnup equations with time-dependent polynomial external feed rates and is much more efficient in improving the accuracy by using substeps, which demonstrates that PRAM is the attractive method for burnup calculations.

Optimizing Deep Learning RNN Topologies on Intel Architecture

Recurrent neural network (RNN) models have been found to be well suited for processing temporal data. In this work, we present an optimized implementation of vanilla RNN cell and its two popular variants: LSTM and GRU for Intel Xeon architecture. Typical implementations of these RNN cells employ one or two large matrix multiplication (GEMM) calls and then apply the element-wise operations (sigmoid/tanh) onto the GEMM results. While this approach is easy to implement by exploiting vendor-optimized GEMM library calls, the data reuse relies on how GEMMs are parallelized and is sub-optimal for GEMM sizes stemming from small minibatch. Also, the element-wise operations are exposed as a bandwidth-bound kernel after the GEMM which is typically a compute-bound kernel. To address this discrepancy, we implemented a parallel blocked matrix GEMM in order to (a) achieve load balance, (b) maximize weight matrix reuse, (c) fuse the element-wise operations after partial GEMM blocks are computed and while they are hot in cache. Additionally, we bring the time step loop in our cell to further increase the weight reuse and amortize the overhead to transform the weights into blocked layout. The results show that our implementation is generally faster than Intel MKL-DNN library implementations, e.g. for RNN, forward pass is up to ~3× faster whereas the backward/weight update pass is up to ~5× faster. Furthermore, we investigate high-performance implementations of sigmoid and tanh activation functions that achieve various levels of accuracy. These implementations rely on minimax polynomial approximations, rational polynomials, Taylor expansions and exponential approximation techniques. Our vectorized implementations can be flexibly integrated into deep learning computations with different accuracy requirements without compromising performance; in fact, these are able to outperform vectorized and reduced accuracy vendor-optimized (Intel SVML) libraries by 1.6–2.6× while speep up over GNU libm is close to two orders of magnitude. All our experiments are conducted on Intel’s latest CascadeLake architecture.

Logistic regression over encrypted data from fully homomorphic encryption

BackgroundOne of the tasks in the 2017 iDASH secure genome analysis competition was to enable training of logistic regression models over encrypted genomic data. More precisely, given a list of approximately 1500 patient records, each with 18 binary features containing information on specific mutations, the idea was for the data holder to encrypt the records using homomorphic encryption, and send them to an untrusted cloud for storage. The cloud could then homomorphically apply a training algorithm on the encrypted data to obtain an encrypted logistic regression model, which can be sent to the data holder for decryption. In this way, the data holder could successfully outsource the training process without revealing either her sensitive data, or the trained model, to the cloud.MethodsOur solution to this problem has several novelties: we use a multi-bit plaintext space in fully homomorphic encryption together with fixed point number encoding; we combine bootstrapping in fully homomorphic encryption with a scaling operation in fixed point arithmetic; we use a minimax polynomial approximation to the sigmoid function and the 1-bit gradient descent method to reduce the plaintext growth in the training process.ResultsOur algorithm for training over encrypted data takes 0.4–3.2 hours per iteration of gradient descent.ConclusionsWe demonstrate the feasibility but high computational cost of training over encrypted data. On the other hand, our method can guarantee the highest level of data privacy in critical applications.

Research on MMPA Algorithm for Solving Detailed Depletion Chain With the NUIT Code

In depletion and transmutation calculation, it is important to solve detailed burnup chains with high computational accuracy and efficiency. This requires the good performance of the burnup algorithms. Nuclide inventory tool (NUIT) is a newly developed nuclide inventory calculation code, which is capable of handling detailed depletion chains by implementing various advanced algorithms. Based on the NUIT code, this paper investigates the accuracy and efficiency of the mini-max polynomial approximation (MMPA) method, and compares it with other burnup solvers in NUIT code. It is concluded that the MMPA method is numerically accurate and efficient for dealing with detailed depletion chains with extremely short half-lived nuclides.

A Robust and Scalable Implementation of the Parks-McClellan Algorithm for Designing FIR Filters

With a long history dating back to the beginning of the 1970s, the Parks-McClellan algorithm is probably the most well known approach for designing finite impulse response filters. Despite being a standard routine in many signal processing packages, it is possible to find practical design specifications where existing codes fail to work. Our goal is twofold. We first examine and present solutions for the practical difficulties related to weighted minimax polynomial approximation problems on multi-interval domains (i.e., the general setting under which the Parks-McClellan algorithm operates). Using these ideas, we then describe a robust implementation of this algorithm. It routinely outperforms existing minimax filter design routines.

Precise and fast computation of generalized Fermi–Dirac integral by parameter polynomial approximation

The generalized Fermi–Dirac integral, Fk(η, β), is approximated by a group of polynomials of β as Fk(η,β)≈∑j=0JgjβjFk+j(η) where J=1(1)10. Here Fk(η) is the Fermi-Dirac integral of order k while gj are the numerical coefficients of the single and double precision minimax polynomial approximations of the generalization factor as 1+x/2≈∑j=0Jgjxj. If β is not so large, an appropriate combination of these approximations computes Fk(η, β) precisely when η is too small to apply the optimally-truncated Sommerfeld expansion (Fukushima, 2014 [15]). For example, a degree 8 single precision polynomial approximation guarantees the 24 bit accuracy of Fk(η, β) of the orders, k=−1/2(1)5/2, when −∞<η≤8.92 and β ≤ 0.2113. Also, a degree 7 double precision polynomial approximation assures the 15 digit accuracy of Fk(η, β) of the same orders when −∞<η≤29.33 and 0≤β≤3.999×10−3. Thanks to the piecewise minimax rational approximations of Fk(η) (Fukushima, 2015 [18]), the averaged CPU time of the new method is roughly the same as that of single evaluation of the integrand of Fk(η, β). Since most of Fk(η) are commonly used in the approximation of Fk(η, β) of multiple contiguous orders, the simultaneous computation of Fk(η, β) of these orders is further accelerated by the factor 2–4. As a result, the new method runs 70–450 times faster than the direct numerical integration in practical applications requiring Fk(η, β).

A Minimax Polynomial Approximation Objective Function Approach for Optimal Design of Power System Stabilizer by Embedding Particle Swarm Optimization

The paper presents a novel approach based on Minimax approximation and evolutionary tool Particle Swarm Optimization (PSO) to fabricate the parameters of Power System Stabilizers (PSSs) for multi machine power systems. The proposed approach employs PSO algorithm for find the setting of PSS parameters. The worth mentioning feature of this work is the formulation of objective function with the help of swing curves interpolation. A novel transformation known as minimax approximation is used for converting the objective into the polynomials of degree one, two and three. To construct the objective function based on interpolation second order sensitivity analysis is performed. The performance of the PSSs is tested under different topological changes, operating conditions and system configurations. Nonlinear simulation reveals that proposed PSSs are effectively deal with local and interarea modes of oscillations. PSS design obtained through lower order polynomial expression of objective function is able to deal with the oscillatory modes efficiently. DOI: http://dx.doi.org/10.11591/telkomnika.v14i2.7602

Numerical solution of matrix exponential in burn-up equation using mini-max polynomial approximation

Nuclear fuel burn-up depletion calculations are essential to compute the nuclear fuel composition transition. In the burn-up calculations, the matrix exponential method has been widely used. In the present paper, we propose a new numerical solution of the matrix exponential, a Mini-Max Polynomial Approximation (MMPA) method. This method is numerically stable for burn-up matrices with extremely short half-lived nuclides as the Chebyshev Rational Approximation Method (CRAM), and it has several advantages over CRAM. We also propose a multi-step calculation, a computational time reduction scheme of the MMPA method, which can perform simultaneously burn-up calculations with several time periods. The applicability of these methods has been theoretically and numerically proved for general burn-up matrices. The numerical verification has been performed, and it has been shown that these methods have high precision equivalent to CRAM.

Barycentric-Remez algorithms for best polynomial approximation in the chebfun system

The Remez algorithm, 75 years old, is a famous method for computing minimax polynomial approximations. Most implementations of this algorithm date to an era when tractable degrees were in the dozens, whereas today, degrees of hundreds or thousands are not a problem. We present a 21st-century update of the Remez ideas in the context of the chebfun software system, which carries out numerical computing with functions rather than numbers. A crucial feature of the new method is its use of chebfun global rootfinding to locate extrema at each iterative step, based on a recursive algorithm combining ideas of Specht, Good, Boyd, and Battles. Another important feature is the use of the barycentric interpolation formula to represent the trial polynomials, which points the way to generalizations for rational approximations. We comment on available software for minimax approximation and its scientific context, arguing that its greatest importance these days is probably for fundamental studies rather than applications.

Robust Solvers for Inverse Imaging Problems Using Dense Single-Precision Hardware

We present an iterative framework for robustly solving large inverse problems arising in imaging using only single-precision (or other reduced-precision) arithmetic, which allows the use of high-density processors (e.g. Cell BE and Graphics Processing Units). Robustness here means linear-convergence even for large problems (billions of variables), with high levels of noise (signal-to-noise levels less than unity). This framework handles problems formulated as quadratic and general non-linear minimization problems. Sparse and dense problems can be treated, as long as there are efficient parallelizable matrix-vector products for the transformations involved. Outer iterations correspond to approximate solutions of a quadratic minimization problem, using a single Newton step. Inner iterations correspond to the estimation of the step via truncated Neumann series or minimax polynomial approximations built from operator splittings. Given the convergence analysis, this approach can also be used in embedded environments with fixed computation budgets, or certification requirements, like real-time medical imaging. We describe a benchmark problem from MRI, and a series of penalty functions suited to this framework. An important family of such penalties is motivated by both bilateral filtering and total variation, and we show how they can be optimized using linear programming. We also discuss penalties designed to segment images, and use different types of a priori knowledge, and show numerically that the different penalties are effective when used in combination.

A Tau method with perturbed boundary conditions for certain ordinary differential equations

Ortiz' recursive formulation of the Lanczos Tau method (TM) is a powerful and efficient technique for producing polynomial approximations for initial or boundary value problems. The method consists in obtaining a polynomial which satisfies (i) aperturbed version of the given differential equation, and (ii) the imposed supplementary conditionsexactly. This paper introduces a new form of the TM, (denoted by PTM), for a restricted class of differential equations, in which the differential equations as well as the supplementary conditions areperturbed simultaneously. PTM is compared to the classical TM from the point of view of their errors: it is found that the PTM error is smaller and more oscillatory than that of the TM; we further find that approximations nearly as accurate as minimax polynomial approximations can be constructed by means of the PTM. Detailed formulae are derived for the polynomial approximations in TM and PTM, based on Canonical Polynomials. Moreover, various limiting properties of Tau coefficients are established and it is shown that the perturbation in PTM behaves asymptotically proprtional to a Chebyshev polynomial.

High-speed double-precision computation of reciprocal, division, square root, and inverse square root

A new method for the high-speed computation of double-precision floating-point reciprocal, division, square root, and inverse square root operations is presented in this paper. This method employs a second-degree minimax polynomial approximation to obtain an accurate initial estimate of the reciprocal and the inverse square root values, and then performs a modified Goldschmidt iteration. The high accuracy of the initial approximation allows us to obtain double-precision results by computing a single Goldschmidt iteration, significantly reducing the latency of the algorithm. Two unfolded architectures are proposed: the first one computing only reciprocal and division operations, and the second one also including the computation of square root and inverse square root. The execution times and area costs for both architectures are estimated, and a comparison with other multiplicative-based methods is presented. The results of this comparison show the achievement of a lower latency than these methods, with similar hardware requirements.

The convergence rate of the minimal residual method for the Stokes problem

Discretisation of the classical Stokes problem gives rise to symmetric indefinite matrices with eigenvalues which, in a precise way, are not symmetric about the origin, but which do depend on a mesh size parameter. Convergence estimates for the Conjugate Residual or Minimum Residual iterative solution of such systems are given by best minimax polynomial approximations on an inclusion set for the eigenvalues. In this paper, an analytic convergence estimate for such problems is given in terms of an asymptotically small mesh size parameter.

Error bounds for a near-minimax approximation

It is known that a near minimax polynomial approximation tof ∈C [−1, 1] is provided by a finite carrier projectionM n fromC[−1, 1] onto the subspace of all polynomials of degree ≤n, such thatM nf is a weighted least squares approximation off on the set consisting of the extreme points of the Chebyshev polynomialT 2n + 1. In this paper, upper bounds for the error ∥f−M n f∥∞ are given in terms of divided differences.

A fast algorithm for linear complex Chebyshev approximations

We propose a new algorithm for finding best minimax polynomial approximations in the complex plane. The algorithm is the first satisfactory generalization of the well-known Remez algorithm for real approximations. Among all available algorithms, ours is the only quadratically convergent one. Numerical examples are presented to illustrate rapid convergence.

Error of truncated Chebyshev series and other near minimax polynomial approximations

It is well known that a near minimax polynomial approximation p is obtained by truncating the Chebyshev series of a function ƒ; after n + 1 terms. It is shown that if ƒ; ϵ C (n + 1)[−1, 1], then ∥ƒ; − p ∥ may be expressed in terms of ƒ; (n + 1) in the same manner as the error of minimax approximation. The result is extended to other types of near minimax approximation.

Polynomial approximation using equioscillation on the extreme points of Chebyshev polynomials

One method of obtaining near minimax polynomial approximation to f ∈ C ( n + 1) [−1, 1] is to choose p ∈ P n such that f − p equioscillates on the point set consisting of the extrema of T n + 1 . It is shown that ∥ f − p∥ may be expressed in terms of f ( n + 1) in the same manner as E n ( f), the error of minimax approximation. The Lebesgue constants are also investigated.

On the foundations of best approximation theory

In this paper we excavate the foundations of best-approximation theory with the tools of Bishop's constructive analysis. We prove a general theorem on existence (computability) of best approximations from a given finite-dimensional linear subspace of a normed space E, and illustrate this with the case where E is uniformly convex. The second part of the paper deals with the characterisation and existence of minimax polynomial approximations to elements of C[0, 1], and with the pointwise continuity of the minimax approximation mapping on this space. In particular, the main application of our general existence theorem answers affirmatively the long-open question: Is there a constructive proof of the existence of minimax polynomial approximations?