The Construction of High Order Convergent Look-Ahead Finite Difference Formulas for Zhang Neural Networks
aa r X i v : . [ m a t h . NA ] A p r The Construction of High Order Convergent Look-AheadFinite Difference Formulas for Zhang Neural Networks
Frank UhligDepartment of Mathematics and Statistics,Auburn University, AL 36849-5310, USA ( [email protected] ) Abstract :
Zhang Neural Networks rely on convergent 1-step ahead finite difference formulas of which veryfew are known. Those which are known have been constructed in ad-hoc ways and suffer from low truncationerror orders. This paper develops a constructive method to find convergent look-ahead finite difference schemesof higher truncation error orders. The method consists of seeding the free variables of a linear system comprisedof Taylor expansion coefficients followed by a minimization algorithm for the maximal magnitude root of theformula’s characteristic polynomial. This helps us find new convergent 1-step ahead finite difference formulasof any truncation error order. Once a polynomial has been found with roots inside the complex unit circle andno repeated roots on it, the associated look-ahead ZNN discretization formula is convergent and can be used forsolving any discretized ZNN based model. Our method recreates and validates the few known convergent formulas,all of which have truncation error orders at most 4. It also creates new convergent 1-step ahead difference formulaswith truncation error orders 5 through 8.
Subject Classifications :
Key Words : finite difference formula, look-ahead difference formula, Taylor expansion, linear systems, freevariable, characteristic polynomial, convergent multistep method, Zhang neural network, truncation error order
Finite difference formulas have a long history of over 200 years in computational mathematics. They came aboutafter the development of Calculus in the late 17th century and were introduced to estimate the behavior and slopeof functions or to approximate areas and volumes. In the differentiation realm, one of the first such formulas thatis still relevant today is Euler’s forward finite difference formula, written here in its symmetric form as ˙ y j ≈ y j +1 − y j − τ (1)where y j = y ( t j ) with t j = t + jτ for a constant step size τ , an initial time instance t and any j ≥ . If wesolve (1) for y j +1 we obtain the symmetric 1-step ahead finite difference formula of Euler y j +1 = y j − + 2 τ ˙ y j . (2)After assembling the y j entries on the left of the equation (2) this leads to y j +1 − y j − = 2 τ ˙ y j . (3)Next we interpret the left hand side of equation (3) as a polynomial p of smallest degree in a variable x where thesubscripts in (3) become the powers of x , namely p ( x ) = x − . p is called the characteristic polynomial of the Frank Uhligdifference equation (3). This process is familiar to anyone who has taken a first course in Numerical Analysis andit is described in every elementary textbook on Numerics. We have explained this fundamental process here in fulldetail because we need to use it repeatedly in the future.The roots of the characteristic polynomial p for Euler’s finite difference formula (2) are +1 and –1. They both lieon the periphery of the unit circle in C and are distinct. Therefore the 1-step Euler method is convergent. Therequirement that all roots of p must lie inside the closed unit disk and no repeated roots may lie on the unit circleis necessary and sufficient for convergence, see e.g. [2, ch. 17.6.2]. Convergent finite difference schemes can beused repeatedly to trace solutions of differential equations in a look-ahead way, albeit for Euler with a low orderof accuracy. This was initially done by Bunse-Gerstner et al in [1] in 1991 for computing time-varying SVDsefficiently via a look-ahead Euler based integrator and subsequently in many other papers.Zhang Neural Networks were first developed early in the new millennium by Yunong Zhang and others, calledZeroing Neural Network then, [18]. The idea was taken up by engineers and implemented in many time-varyingapplications with well over three hundred articles published mostly in engineering and applied math journals. Re-cently there have been impulses from numerical analysts, but ZNN is still largely unknown in numerical circles.ZNN based methods have been used to find time-varying reciprocals, square roots, generalized inverses, pseu-doinverses, and for solving linear, Sylvester, Lyapunov and other equations or inequalities, as well as for matrixeigenvalues and eigenvectors and almost anything matrix related; in optimal design and control and for minimiza-tions and so forth. The list of real-world applications goes on and on, from robotics to autonomous cars andexperimental aircrafts where sensor data arrives at frequencies such as 50 Hz and an objective must be met accu-rately, 1-step ahead in real-time and for real-world situations. ZNN methods can be easily implemented in on-chipdesigns for practical control applications, see e.g. [17]. They currently have many industrial uses and give systemsand machines improved performance. See e.g. [4, 5, 10, 13, 18, 19, 20, 21, 22, 24] for a short glimpse of theirvast potential with time-varying systems. Most recently there have been substantial advances on the numericalbehavior of continuous-time ZNN methods for matrix problems such as by Lin Xiao et al in [16] on the stabilityand robustness in control application, in [15] on robot applications, and on general time-varying matrix inverseproblems when using time-varying decay constants η ( t ) > in [12, 11] for example.Zhang Neural Networks (ZNN) are designed for and most efficiently used to solve time-varying multi-dimensionalequations f ( t ) = b ( t ) predictively with high accuracy and quickly in real-time. There they use look-ahead finitedifference formulas to solve a problem specific error differential equation. All ZNN methods are based on the errorequation e ( t ) = f ( t ) − b ( t ) ! = 0 and the stipulation that e ( t ) should decay exponentially fast to zero, i.e., ˙ e ( t ) = − ηe ( t ) for η > . (4)When associated with a given continuous-time f ( t ) = b ( t ) model, discretized ZNN methods can easily solve theassociated discrete time-varying problem f ( t k ) = b ( t k ) where the time steps t k are equidistant and the input datais derived from repeated sensor readings. The discretized method predicts the model’s solution at time t k +1 from acertain subset of the previous iterates f ( t j ) and b ( t j ) with j ≤ k and does so shortly after time t k by using certainderivatives that the error DE (4) requires. Discrete ZNN methods must construct the next iterate well before thenext time instance t k +1 arrives. They succeed with high accuracy and a truncation error order of O ( τ m +1 ) if thechosen 1-step ahead discretization formula has truncation error order m + 1 for the constant sampling gap τ .Section 2 below describes how to set up the mechanics for finding high order convergent look-ahead finitedifference formulas via Linear Algebra and elementary Matrix Theory. Section 3 then describes a characteristicpolynomial root minimization process that helps us find look-ahead finite difference schemes which satisfy theconvergence conditions. Section 4 provides a list of new and high truncation order convergent 1-step ahead finitedifference schemes of truncation error orders up to O ( τ ) , as well as open problems. Several of our newly foundconvergent 1-step ahead finite difference formulas are tested regarding their accuracy and convergence behavior in[8] on the parameter-varying complex matrix field of values problem.onvergent look-ahead Difference Formulas 3 A recent literature search for known convergent look-ahead finite difference formulas found less than half a dozensuch formulas, all of which had truncation error orders less than or equal to 4. This is so despite hundreds ofpotentially known look-ahead discretization formulas in the literature. Most of them fail the characteristic rootscondition and therefore they are unusable for ZNN type methods, such as the look-ahead but unstable formulas(27) to (30) in [4]. With repeated roots on the unit circle, oscillations will eventually set in; with roots outside theunit circle, divergence to infinity will occur whenever corresponding fundamental solutions creep into the currentstates. With this paper we more than double the range of available convergent look-ahead methods from truncationerror orders 2, 3, and 4, up to error orders 5, 6, 7 and 8. The new high error order formulas speed up convergenceand improve accuracy to near machine constant error levels when properly implemented.This linear algebraic section develops the basis for an algorithm to determine look-ahead finite difference formulasof any truncation order, regardless of convergence or not. The computed look-ahead methods will then be used inthe next section to act as initial guesses or seeds for starting a roots minimization process that may find convergentlook-ahead schemes or it might not. Note that the previously found convergent methods were all (except for Euler)computed by ad hoc methods with lucky guesses and clever schemes. Here the process is formalized and developedinto a computer code that needs no luck, no sweat and no tears.Let us consider a discrete time-varying state vector x j = x ( t j ) = x ( j · τ ) for a constant sampling gap τ and j = 0 , , , ... and write out ℓ + 1 explicit Taylor expansions for x j +1 , x j − , ..., x j − ℓ around x j as follows: x j +1 = x j + τ ˙ x j z }| { + τ
2! ¨ x j + τ ... x j ... + τ m m ! m ˙ x j + O ( τ m +1 ) (5) x j − = x j − τ ˙ x j + τ
2! ¨ x j − τ ... x j ... + ( − m τ m m ! m ˙ x j + O ( τ m +1 ) (6) x j − = x j − τ ˙ x j + (2 τ )
2! ¨ x j − (2 τ ) ... x j ... + ( − m (2 τ ) m m ! m ˙ x j + O ( τ m +1 ) (7) x j − = x j − τ ˙ x j + (3 τ )
2! ¨ x j − (3 τ ) ... x j ... + ( − m (3 τ ) m m ! m ˙ x j + O ( τ m +1 ) (8)... (9) x j − ℓ = x j − ℓτ ˙ x j + ( ℓτ )
2! ¨ x j − ( ℓτ ) ... x j ... + ( − m ( ℓτ ) m m ! m ˙ x j | {z } + O ( τ m +1 ) (10)Each right hand side of the Taylor expansion rows or equations above contains m + 2 terms. The central under-and overbraced m − ’column terms’ on the right hand side of the equal signs each contain a factor of identicalpowers of τ and identical varying order partial derivatives of x j . Our interest lies only in the remaining ’rationalnumber’ factors in the third through ( m + 1) st ’columns’ on the right hand side of equations (5) through (10),i.e., for the moment we omit the powers of τ and the derivatives r ˙ x j for r = 2 , ..., m throughout what immediatelyfollows. If we can find a linear combination of the ℓ + 1 equations (5) through (10) that makes the ’braced’ m − number terms in each of these ’columns’ disappear or become 0, we have found an equation for x j +1 in terms ofthe already known values of x j , x j − , ..., x j − ℓ , the derivative ˙ x j with an overall error term of order O ( τ m +1 ) .These are the only items that are left once the braced region’s linear row combination has become 0. To zero out Frank Uhligthe braced region, we now collect the relevant rational numbers factors in the constant matrix A ℓ +1 ,m − A =
12! 13! 14! · · · m !12! −
13! 14! · · · ( − m m !2 −
3! 2 · · · ( − m m m ! ... ... ... ℓ − ℓ ℓ · · · ( − m ℓ m m ! ∈ R ℓ +1 ,m − . (11) A ’s entry a u,v in row u and column v is ( − v +1 ( u − v +1 ( v + 1)! for ≤ u ≤ ℓ + 1 and a ,v = 1( v + 1)! for all v = 1 , ..., m − . The complete over- and underbraced summed terms in equations (5) through (10) has the matrixtimes vector product form A ℓ +1 ,m − · taudx =
12! 13! 14! · · · m !12! −
13! 14! · · · ( − m m !2 −
3! 2 · · · ( − m m m ! ... ... ... ℓ − ℓ ℓ · · · ( − m ℓ m m ! τ ¨ x j τ ... x j τ ˙ x j ... τ m − m − ˙ x j τ m m ˙ x j (12)where the vector taudx ∈ R m − contains the increasing powers of τ multiplied by the respective higher deriva-tives of x j as entries. If we can find a left kernel row vector x ∈ R ℓ +1 for A with x · A = o m − , the zerorow vector in R m − , then x · A · taudx = 0 ∈ R as well. If we then form the linear combination of the ℓ + 1 equations in rows (5) through (10) as prescribed by the coefficients of x , the linear combination of the sums inthe under- and overbraced columns in (5) through (10) will vanish and we obtain a single look-ahead formula, in-volving only x j +1 , x j , x j − , ..., x j − ℓ , and ˙ x j with an error term of order O ( τ m +1 ) . A non-zero left kernel vector x for A ℓ +1 ,m − exists as soon as the number of rows ℓ +1 of A exceeds the number m − of columns of A ℓ +1 ,m − .A left null row vector x for A is a right null column vector for A T when transposed and vice versa. Columnnull vectors can be found from a reduced row echelon form reduction R of A T easily by setting the free variablesof R m − ,ℓ +1 equal to a nonzero seed vector y of length ℓ + 1 − ( m −
1) = ℓ − m + 2 . In this case – which willbe called the regular case from now on, R is a single row block matrix with the identity matrix I m − appearingin the first block position because the coefficient matrix A always has full rank m − and the linear system isunderdetermined. And generally, a dense m − by m − matrix B appears in the second block position of R , i.e., R has the form R = (cid:0) I m − , B m − ,m − (cid:1) m − , m − . (13)This dimensional situation works very well here. In fact the method works for any matrix A k + s,k and anynonzero seed vector y ∈ R s . Any nonzero seed vector y ∈ R s spawns a null vector q ∈ R k + s for A T as q = [-R*[zeros(k,1);y];y] in Matlab notation. We then replace the vector q by q/q (1) in order to arriveat a normalized characteristic polynomial p for the associated convergent look-ahead finite difference equation.Once such a null vector q ∈ R m − has been computed from a seed vector y ∈ R m − in the regular A m − ,m − case where ℓ +1 = 2( m − , our algorithm forms the specific linear combination of the set of equations (5) through(10) that q suggests in order to zero out all contributions from the entries in the under- and overbraced third through m + 1 st columns on the right hand side of equations (5) to (10). Then we separate the 1-step ahead state x j +1 onthe left hand side with likewise accumulations for the current and earlier states x j , x j − , ..., x j − m +3 and the firstonvergent look-ahead Difference Formulas 5time derivative ˙ x j according to q on the right hand side of the linearly combined equations. This process yields afinite difference multistep formula for m − equidistant state vectors. Its characteristic polynomial is the nor-malized polynomial p = [1;-sum(q);q(2:2(m-1))] ∈ R m − , again in Matlab notation. The polynomial p describes a 1-step-ahead difference equation for m − contiguous instances, i.e., a (2 m − -IFD formula in ourabbreviation to denote ’ (2 m − -Instance Difference Formulas’ . x j +1 + p x j + ... + p m − x j − m +3 = c ˙ x j . (14)We solve (14) for x j +1 by rearranging terms and obtain the associated look-ahead rule x j +1 = − ( p x j + ... + p m − x j − m +3 ) + c ˙ x ˜ j (15)where we have incorporated the linear combination of the first derivative ˙ x j terms in equations (5) through (10) inthe constant c = p . ∗ [1; 0; − (1 : 2 m − ′ ] (in Matlab notation). Here p (with leading coefficient p normalizedto 1) is in row vector form to create the dot product c . And ˙ x ˜ j is a backward approximation of the derivative ˙ x j attime t j of sufficiently high truncation error order that can by computed from previous x .. state data. To use formula(15) in the A m − ,m − regular case we clearly need m − starting values for the x .. terms. And we need k + s starting values x .. in the more general A k + s,k situation.Our next task is to find convergent finite difference 1-step ahead formulas of a form such as (15). Convergence ofsuch multistep formulas depends on the lay of their characteristic polynomial’s roots: these must all lie inside theunit circle in C and if there are roots on the unit circle, those must be simple. See e.g. [2, Ch. 17.6.2 and Definition17.17, p. 475]. How to achieve convergence behavior here is the subject of our next Section. In Section 2 we have described how to start from a short seed vector and construct a 1-step ahead finite differencescheme for possible use in general ZNN processes. The behavior of the associated characteristic polynomial deter-mines convergence or divergence. If the constructed scheme is not convergent for the chosen seed y then there areroots of its characteristic polynomial that exceed 1 in magnitude or there are repeated roots on the unit circle. Inour experience, the former, i.e., roots outside the unit circle, is seemingly always the cause for non-convergence;we have never encountered the latter. Note that the matrix A k + s,k has very special rational number entries ofinteger powers divided by factorials. Therefore the reduced row echelon form R of A T contains mostly integers orrational numbers with small integer denominators in its second block B as defined in (13).In this section we propose a minimizing algorithm for the maximal magnitude characteristic polynomial root interms of a given seed vector y . For A k + s,k the seed vector space is R s and any seed y therein spawns a uniquelook-ahead finite difference scheme and an associated characteristic polynomial as was shown in Section 2. How-ever, the set of possible characteristic polynomials itself is not a linear space since sums of such polynomials mayor may not be representatives of look-ahead difference schemes. Therefore we can only vary the seed and not theintermediate polynomials in the minimization process and we will have to search indirectly for a characteristicpolynomial with smaller maximal magnitude root in a neighborhood of the starting seed y . We implement thisminimization process by using the multidimensional built-in Matlab fminsearch minimizer function until wehave found a seed with an associated characteristic polynomial that is convergent or there is no convergent suchformula from the chosen seed. fminsearch uses the Nelder-Mead downhill simplex method [6] that finds localminima for non-linear functions such as ours without using derivatives. It mimics the method of steepest descendand performs a local minimizing search via multiple function evaluations. The main task is to discover generatingseed vectors y that can start the minimizing iteration and find a characteristic polynomial that is convergent ac-cording to the convergent root stipulations. Our seed selection process is currently based on random entry seedssince we know of no better way.For the general A k + s,k case and after many different approaches for choosing our random entry seed vectors y ∈ R s , we decided to start from seeds y with normally distributed entries and run fminsearch to try and find Frank Uhliga local minimum of the maximal magnitude root of the associated characteristic polynomial. Our minimizationalgorithm runs through a double do loop. An outside loop for a number (5 to 20 or 100 ...) of random entry startingseed vectors and an inner loop for a number (4 to 7 or 15 ...) of randomized restarts from a previously computed fminsearch polynomial that is non-convergent. We project its seed onto a point with newly randomized entriesnearby and use this new seed for a subsequent inner loop run several (2, 5 or 8) times.The whole MATLAB code fits onto 80 lines of code, plus 30 lines of comments and two dozen lines of knownconvergent 1-step ahead discretization formulas of varying truncation error orders between O ( τ ) to O ( τ ) . Pre-viously convergent look-ahead methods were completely unknown for truncation error orders above O ( τ ) . Nowwe can compute many convergent polynomials of higher error orders quickly. Here is the complete short list of the six known 1-step ahead discretization formulas and their properties orderedby ascending truncation error order. (A) Symmetric Euler Differentiation and Discretization Formula with ZNN truncation error order O ( τ ) : ˙ y j = 12 τ y j +1 − τ y j − + O ( τ ) or y j +1 = y j − + O ( τ ) + ... problem specific terms from model’s right hand side and ˙ x j Characteristic Polynomial : p ( x ) = x − Formally and according to multistep theory, the Euler formula should lead to a convergent look-ahead discretiza-tion formula, but Euler is not convergent as such in practice. Why so is a surprising mystery. (B) 4-IFD Formula from [5, equations (10), (12)] with ZNN truncation error order O ( τ ) : (10) ˙ y j = 1 τ y j +1 − τ y j + 1 τ y j − − τ y j − + O ( τ ) or (12) y j +1 = 32 y j − y j − + 12 y j − + O ( τ ) + ... problem specific terms from model’s rhs and ˙ x j Characteristic Polynomial : p ( x ) = 2 x − x + 2 x − , (not normalized) (C) 4-IFD Formula from [5, equation (11)] with ZNN truncation error order O ( τ ) : (11) ˙ y j = 35 τ y j +1 − τ y j − τ y j − − τ y j − + O ( τ ) or y j +1 = 12 y j + 13 y j − + 16 y j − + O ( τ ) + ... problem specific terms from ...Characteristic Polynomial : p ( x ) = 6 x − x − x − , (not normalized) (D) FIFD Formula from [4, equations (14), (21)] with ZNN truncation error order O ( τ ) : (14) ˙ y j = 58 τ y j +1 − τ y j − τ y j − − τ y j − + O ( τ ) or (21) y j +1 = 35 y j + 15 y j − + 15 y j − + O ( τ ) + ... problem specific terms from ...Characteristic Polynomial : p ( x ) = 5 x − x − x − , (not normalized) (E) 5-IFD Formula from [5, equations (23), (27)] with ZNN truncation error order O ( τ ) : (23) ˙ y j = 49 τ y j +1 + 118 τ y j − τ y j − − τ y j − + 19 τ y j − + O ( τ ) or (27) y j +1 = − y j + 34 y j − + 58 y j − − y j − + O ( τ ) + ... problem specific terms from ...Characteristic Polynomial : p ( x ) = 8 x + x − x − x + 2 , (not normalized)onvergent look-ahead Difference Formulas 7 (F) 6N τ CD Formula from [9, equations (16), (18)] with ZNN truncation error order O ( τ ) : (16) ˙ y j = 1324 τ y j +1 − τ y j − τ y j − − τ y j − − τ y j − + 112 τ y j − + O ( τ ) or (18) y j +1 = 613 y j + 213 y j − + 413 y j − + 313 y j − − y j − + O ( τ ) + ... problem specific terms ...Characteristic Polynomial : p ( x ) = 13 x − x − x − x − x + 2 , (not normalized)Next we mention a list of new convergent look-ahead discretisation formulas that is detailed within our MAT-LAB codes. We start with some results obtained using our runconv1step.m file [7] in Matlab. There are 4integer inputs to runconv1step.m : the first input indicates how many outer loop runs are desired and the sec-ond input indicates how many separate repeats with altered seed inputs should be performed in an inner do loop foreach outer run. A call of runconv1step(40,10,k,s) would thus require 40 outer loop runs with 10 separateseeds each (400 runs of fminsearch in total). This call tries to find convergent polynomials of truncation errororders k + 2 , i.e., polynomials of degree k + s with roots properly inside the closed unit disk in C . Here k can beany integer less than 6 and s should be at least equal to k so that the rational entry matrix A k + s,k has a nontrivialleft nullspace. Here s is the number of real entries in the seed vector y and our highest tested-for truncation errororder O ( τ k +2 ) is k + 2 = 6 + 2 = 8 for k = 6 . We see no real need to go beyond k = 6 since a truncation errororder of (50 Hz ) = 0 . ≈ . · − seems close enough to the machine constant to wonder about furtherimprovements.The MATLAB output for validating the known convergent 5-IFD formula (E) above from [5] which is listedin the examples list in our code runconv1step.m [7] is as follows with y denoting the seed vector of length s that was used. format short, tic, TPOLY = runconv1step(1,1,2,2), tocTruncation_error_order =4y = -5 2TPOLY = 1.0000 0.1250 -0.7500 -0.6250 0.2500 0 -0.0000 0.9025 2.2500Elapsed time is 0.008423 seconds. The first k + s + 1 entries in TPOLY are the computed coefficients of the normalized convergent polynomial p ( x ) = x + 0 . x − . x − . x + 0 . of degree k + s = 4 in decreasing exponent order. These arefollowed by a data separating zero, the deviation of the maximal magnitude root of p from 1 (which is nearly zero),the magnitude of the second largest magnitude root of p and finally the coefficient of τ that is to be used inside theZNN discretization. The reader can sample any of the ZNN papers in the bibliography to learn how to set up andimplement a discrete ZNN method from any known convergent polynomial.Next a similar example for a convergent truncation error order 5 formula from the seed y = [a,110,-40] withvariable constant ≤ a ≤ . , also given in our formula list inside runconv1step.m . We have set a = 1 toobtain the output below. format rat, tic, TPOLY = runconv1step(1,1,3,3), tocTruncation_error_order =5y = 1 110 -40TPOLY =1 80/237 -182/237 -206/237 1/237 110/237 -40/237 0 0.0000 446/465 196/79Elapsed time is 0.007906 seconds. Frank UhligNote that the resulting normalized characteristic polynomial has again only rational coefficients. This occurrenceis rather rare. If we run runconv1step(20,10,5,5) with ·
10 = 200 individual searches we typicallycapture 4 to 6 convergent polynomials with discrete ZNN model truncation error orders O ( τ ) . When doing thesame with k = 4 and s = 4 in 200 searches, our code discovers around 40 convergent polynomials with asso-ciated truncation error order O ( τ ) in ZNN applications. It seems advantageous to increase s beyond s = k tofind more higher truncation error order polynomials quickly. For example runconv1step(20,10,5,6) with s = k + 1 = 6 exhibits around a dozen good polynomials and runconv1step(20,10,5,7) with s set to k +2 = 7 around 18. These success numbers are influenced by the random nature of our seeds y of any fixed length s ≥ k . Note that if the call of runconv1step(runs,jend,k,s) generates feasible polynomials, these willall be of degree k + s . Therefore an implementation inside any discrete ZNN model must compute k + s startingvalues before an associated 1-step forward discrete ZNN iteration can be run.We have found very few high order convergent polynomials with integer coefficients. The usual output TPOLY from runconv1step(...,...,k,s) comes in 16 digit exponent 10 Matlab notation. The look-ahead con-vergent finite difference formulas of high truncation error orders O ( τ ) to O ( τ ) at the bottom of our our list in runconv1step(...,...,k,s) can be generated from the given seed information there. We have not testedthe quality of the newly computed convergent polynomials and whether there are any noticeable differences inreal-world problems between them. For example we do not know whether those polynomials with well separatedlargest and second largest magnitude roots perform better than those with multiple very near 1 magnitude roots. Tostudy this issue, we have included the magnitude of the second largest root in the last but one column of TPOLY .This paper describes how we can now create high error order convergent look-ahead finite difference schemes forany k and s ≥ k . Our newly found k_s designated polynomials , , and of truncation orders 5through 7 were tested recently on the matrix field of values problem and compared in [8].We conclude with a more theoretical aspect and an open question encountered with finding convergent poly-nomials for 1-step ahead ZNN processes. Remark and
Open Question :
Every polynomial that we have constructed from any seed vector y ∈ R s by our method has had at least oneroot on the unit circle (within − numerical accuracy). This is so even for non-convergent polynomials withsome roots outside the unit disk. Is it true in general that all such Taylor expansion matrix A based polynomialshave at least one root on the unit circle in C ? Or are there some such convergent polynomials whose roots alllie inside the open unit disk. Our method to construct convergent look-ahead finite difference equations hinges on two separate ideas and subse-quently two constructive steps. The first problem, part one, is essentially linear. We want to eliminate the second to m th derivatives in the Taylor expansion equations (5) to (10) around x j and find a certain linear combination of the ℓ + 1 equations that can give us a candidate finite difference scheme. This linear first branch of our task starts froma short seed vector y . It completes the seed of length s to a full difference equation of length k + s and the associ-ated characteristic polynomial coefficients. The resulting difference equation relates the 1-step ahead x j +1 value toearlier x k values with k ≤ j and to the derivative at x j with an error term of order O ( τ m +1 ) . The second problem,part two, is highly nonlinear. It tries to select those finite difference equations whose characteristic polynomialssatisfy certain root conditions for convergence. To solve the second, the nonlinear polynomial roots problem, wehave chosen a multidimensional minimization function that starts from a candidate finite difference scheme an itsassociated set of coefficients and varies the original seed with an eye on minimizing the largest magnitude rootof the associated varying characteristic polynomial. It is the nature of the very first seed from the linear part onethat eventually determines the convergence qualities of the method at the end of each minimizing part two process.Only then do we know.In most of our trial test runs from a seed to a possibly convergent finite difference scheme, the original maximalmagnitude characteristic polynomial root does not dip down to 1 or below. Instead the maximal magnitude rootusually settles around 1.008.., 1.01.. or larger and does no longer budge, indicating that the original seed vector isonvergent look-ahead Difference Formulas 9not allowing our algorithm to find a usable difference scheme and that we must end this run. Clearly not all seedscan solve our problem, their surrounding maximal root ’valleys’ simply may not dip low enough to be of use to us.But often enough, our randomized seed minimization algorithm leads to success in finding convergent look-aheadfinite difference formulas with truncation error orders up to 8 where none were previously known. References [1] A. Bunse-Gerstner, R. Byers, V. Mehrmann, N. K. Nichols, Numerical computation of an analytic singularvalue decomposition of a matrix valued function, Numerische Mathematik, 60 (1991), 1 - 39.[2] G. Engeln-M¨ullges, F. Uhlig,
Numerical Algorithms with C , Springer, 1996.[3] J. C. Lagarias, J. A. Reeds, M. H. Wright, P. E. Wright, Convergence properties of the Nelder-Mead simplexmethod in low dimensions, SIAM J of Optimization, 9 (1998), 112 - 147.[4] J. Li, M. Mao, F. Uhlig, Y. Zhang, Z-type neural-dynamics for time-varying nonlinear optimization undera linear equality constraint with robot application, Journal of Computational and Applied Mathematics, 327(2018), 155 - 166.[5] J. Li, M. Mao, F. Uhlig, Y. Zhang, A 5-instant finite difference formula to find discrete time-varying general-ized matrix inverses and scalar reciprocals, 17 p., submitted[6] J. A. Nelder, R. Mead, A simplex method for function minimization, Computer Journal, 7 (1965), 308 - 313.[7] Frank Uhlig, Matlab code runconv1step(runs,jend,k,s) is available at .[8] Frank Uhlig, Zhang neural networks for fast and accurate computations of the field of values, submitted, 13p. Preprint available as FOVviaZNNa.pdf at .[9] Binbin Qiu, Yunong Zhang, Jinjin Guo, Zhi Yang, Xiaodong, New five-step DTZD algorithm for futurenonlinear minimization with quartic steady-state error pattern, Numerical Algorithms, (2018), in print,https://doi.org/10.1007/s11075-018-0581-4.[10] B. Qiu, Y. Zhang, Z. Yang, New Discrete-Time ZNN Models for Least-Squares Solutions of Dynamic Lin-ear Equation Systems With Time-Varying Rank-Deficient Coefficients, IEEE Trans Neural Networks andLearning Systems, to appear (2019), 12 p.[11] X.-Z. Wang, P. S. Stanimirovic, Y. Wei, Complex ZFs for computing time-varying complex outer inverses,Neurocomputing, 275 (2018), 983 - 1001.[12] X.-Z. Wang, Y. Wei, P. S. Stanimirovic, Complex neural network models for time-varying Drazin nverse,Neural Computation, 28 (2016), 2790 - 2824.[13] L. Xiao, B. Liao, A convergence-accelerated Zhang neural network and its solution application to Lyapunovequation, Neurocomputing, 193 (2016), 213 - 218.[14] Lin Xiao, Bolin Liao, Shuai Li, Ke Chen, Nonlinear recurrent neural networks for finite-time solution ofgeneral time-varying linear matrix equations, Neural Networks, 98 (2018), 102 - 113.[15] Lin Xiao, Zhijun Zhang, Zili Zhang, Weibing Li, Shuai Li, Design, verification and robotic application ofa novel recurrent neural network for computing dynamic Sylvester equation, Neural Networks, 105 (2018),185 - 196.0 Frank Uhlig[16] Lin Xiao, Shuai Li, Faa-Jeng Lin, Zhiguo Tan, Ameer Hamza Khan, Zeroing neural dynamics for controldesign: comprehensive analysis on stability, robustness, and convergence speed, IEEE Transactions on In-dustrial Informatics, in print (2019), DOI 10.1109/TII.2018.2867169.[17] Y. Zhang, D Guo,
Zhang Functions and Various Models , Springer, 2015.[18] Y. Zhang, D. Jiang, J. Wang, A recurrent neural network for solving Sylvester equation with time-varyingcoefficients, IEEE Transactions on Neural Networks, 13 (5) (2002), 1053 - 1063.[19] Y. Zhang, Z. Ke, P. Xu, C. Yi, Time-varying square roots finding via Zhang dynamics versus gradient dynam-ics and the former’s link and new explanation to Newton-Raphson iteration, Information Processing Letters,110 (24) (2010), 1103 - 1109.[20] Y. Zhang, F. Li, Y. Yang, Z. Li, Different Zhang functions leading to different Zhang-dynamics modelsillustrated via time-varying reciprocal solving, Applied Mathematical Modeling, 36 (9) (2012), 4502 - 4511.[21] Y. Zhang, M. Yang, J. Li, L. He, S. Wu, ZFD formula 4IgSFDY applied to future minimization, PhysicsLetters A, 381 (19) (2017), 1677 - 1681.[22] Y. Zhang, M. Yang, C. Li, F. Uhlig, H. Hu, New continuous ZD model for computation of time-varyingeigenvalues and corresponding eigenvectors, submitted, 17 p.[23] Y. Zhang, C. Yi,