Solving Nonlinear Algebraic Systems Using Artificial

and one-dimensional eigenvalues coincide with the zerodimensional and one-dimensional roots of the nonlinearsystem under consideration. On the other hand, at the second stage, the spectral problem for D(λ, µ) is solved, i.e.all zero-dimensional eigenvalues of D(λ, µ) as well as aregular polynomial matrix pencil whose spectrum gives allone-dimensional eigenvalues of D(λ, µ) are found. Regarding the second method, it is based to the factorization ofF (λ, µ) into irreducible factors and to the estimation of theroots (µ, λ) one after the other, since the resulting polynomials produced by this factorization are polynomials of onlyone variable.The last family of methods mentioned in this section isthe one proposed by Emiris, Mourrain and Vrahatis [9].These methods are able to count and identify the roots ofa nonlinear algebraic system based to the concept of topological degree and by using bisection techniques.3 ANNs as linear systems solversAn effective tool for solving linear as well as nonlinearalgebraic systems, is the artificial neural networks (ANNs).Before the description of their usage in the nonlinear case,let us briefly describe how they can be used to solve linearsystems of m equations with n unknowns in the form Az B [11]. In this approach the design of the neural networkallows the minimization of the Euclidean distance r Az b or the equivalent residual r Cz d where C AT Aand d AT b. The proposed neural network is a Hopfieldnetwork with a Lyapunov energy function in the formZnnnnXX1X X1 vi 1E Tij vi vj bi vi gi (v)dv (3)2 i 1ri 1i 1 i 0j6 i 1where bi di /Cii is the externally applied stimulation,vi zi gi (ui ) ui is the neuron output, ui is theneuron state, gi (ui ) is the neuron activation functions and½ Cij /Cii i 6 jTij Tji (4)0i j(i, j 1, 2, . . . , n) are the synaptic weights. The neuralnetwork has been designed in such a way, that after training, the outputs of the network to be the roots of the linearsystem under consideration. This result has been theoretically established by proving that the Lyapunov function ofthe proposed neural network is minimized by the roots ofthe given linear system.4 Theory of algebraic nonlinear systemsAccording to the basic principles of the nonlinear algebra [12], a complete nonlinear algebraic system of n polynomial equations with n unknowns z (z1 , z2 , . . . , zn )is identified completely by the number of equations, n andtheir degrees (s1 , s2 , . . . , sn ), it is expressed mathematically asAi ( z) nXj1 j2 .jsiAizj1 zj2 . . . zjsi 0(5)j1 ,j2 ,.,jsi(i 1, 2, . . . , n), and it has one non-vanishing solution(i.e. at least one jj 6 0) if and only if the equationj1 j2 .jsi s1 ,s2 ,.sn {Ai} 0 holds. In this equation, thefunction is called the resultant and it is a straightforwardgeneralization of the determinant of a linear system. The resultant is a polynomial of the coefficients of A of degree¶n µYXds1 ,s2 ,.,sn degA s1 ,s2 ,.,sn sj(6)i 1j6 iWhen all degrees coincide, i.e. s1 s2 · · · sn s, the resultant n s is reduced to a polynomial of degreedn s degA n s nsn 1 and it is described completelyby the values of the parameters n and s.To understand the mathematical description of a complete nonlinear algebraic system let us consider the case ofn s 2. By applying the defining equation, the ith equation (i 1, 2, . . . , n) of the nonlinear system is estimatedasAi ( z) 2 X2XAji 1 j2 zj1 zj2 j1 1 j2 122 221122 A11i z1 (Ai Ai )z1 z2 Ai z2(7)where the analytical expansion has been omitted for thesake of brevity. The same equation for the case n 3 ands 3 gets the form3 X3 X3XAi ( z) Aji 1 j2 j3 zj1 zj2 zj3 j1 1j2 1j3 1121211222 3333 33 A111i z1 Ai z2 Ai z3 (Ai Ai 131311 22 A113 A112i Ai )z1 z3 i )z1 z2 (Ai3312212122 A313 (A221ii Ai Ai )z1 z2 (Ai13323322333232 Ai )z1 z3 (Ai Ai Ai )z2 z3 223322 2321231 (A232i Ai Ai )z2 z3 (Ai Ai 132231123 A312(8)i Ai Ai Ai )z1 z2 z3From the above equations, it is clear that the coefficients ofthe matrix A which is actually a tensor for n 2 are not allindependent each other. More specifically, for the simplecase s1 s2 · · · sn s, the matrix A is symmetricin the last s contravariant indices. It can be proven thatsuch a tensor has only nMn s independent coefficients, withMn s (n s 1) /(n 1)!s!.

Even though the notion of the resultant has been definedfor homogenous nonlinear equations, it can also describenon-homogenous algebraic equations as well. In the generalcase, the resultant , satisfies the nonlinear Craemer rule s1 ,s2 ,.,sn {A(k) (Zk )} 0(9)where Zk is the kth component of the solution of the nonhomogenous system, and A(k) the kth column of the coefficient matrix, A.Since in the next sections neural models for solving thecomplete 2 2 and 3 3 nonlinear algebraic systems areproposed, let us describe their basic features for the simplecase s1 s2 s. The complete 2 2 nonlinear system isdefined as A(x, y) 0 and B(x, y) 0 whereA(x, y) B(x, y) sXαk xk y s k αsk 0sXsY(x λj y) y s Ã(t) (10)j 1βk xk y s k βsk 0sY(x µj y) y s B̃(t) (11)(α11 X 2 α14 X β1 )z 2 (α13 X α15 )yz α12 y 2 0 (18)(α21 X 2 α24 X β2 )z 2 (α23 X α25 )yz α22 y 2 0 (19)andα11 x2 (α13 Y α14 )xz (α12 Y 2 α15 Y β1 )z 2 0 (20)α21 x2 (α23 Y α24 )xz (α22 Y 2 α25 Y β1 )z 2 0 (21)has the form (z, y) (1, Y ) for the first system, and(x, z) (X, 1) for the second system. But this implies,that the corresponding resultants vanish, i.e., the X variablesatisfies the equation α11 X 2 α14 X β1 ; α13 X α15 ; α12 ; 0 0; α11 X 2 α14 X β1 ; α13 X α15 ; α12 α21 X 2 α24 X β2 ; α23 X α25 ; α22 ; 0 0 (22) 0; α21 X 2 α24 X β2 ; α23 X α25 ; α22 j 1with t y/x, x z1 and y z2 . The resultant of thissystem, has the form¶ss µYY11 (αs βs )s (λi µj ) (α0 β0 )s (12)µj λii,j 1i,j 1(for sake of simplicity we used the notation 2 s {A, B}) and it can be expressed as the determinant ofthe 2s 2s matrix of coefficients. In the particular case ofa linear map (s 1), this resultant reduces to the determinant of the 2 2 matrix, and therefore, it has the form 2 1 {A} kα1 ; α0 ; β1 ; β0 k. On the other hand, for s 2(this is a case of interest in this project) the homogenousnonlinear algebraic system has the formα11 x2 α13 xy α12 y 2α21 x2 α23 xy α22 y 2with a resultant in the form α11 0 2 2 2 {A} α21 0To solve this system, we note that if (X, Y ) is the desiredsolution, then, the solution of the homogenous systemsα13α11α23α21 0 0α12α13α22α230α120α22(13)(14) (15)In complete accordance with the theory of the linear algebraic systems, the above system has a solution if the resultant satisfies the equation 0. Regarding the nonhomogenous complete 2 2 nonlinear system, is can bederived from the homogenous one, by adding and the linearterms; it therefore, has the formα11 x2 α13 xy α12 y 2 α14 x α15 y β1α21 x2 α23 xy α22 y 2 α24 x α25 y β2(16)(17)while, the Y variable satisfies the equation α11 ; α13 Y α14 ; α12 Y 2 α15 Y β1 ; 0 0; α11 ; α13 Y α14 ; α12 Y 2 α15 Y β1 α21 ; α23 Y α24 ; α22 Y 2 α25 Y β2 ; 0 0; α21 ; α23 Y α24 ; α22 Y 2 α25 Y β1 0 (23) (in the above equations the symbol ”;” is used as columnseparator in the tabular environment). Therefore, the variables X and Y got separated, and they can be estimatedfrom separate algebraic equations. However, these solutionsare actually correlated: the above equations are of the 4thpower in X and Y respectively, but making a choice of oneof the four X 0 s, one fixes associated choice of Y . Thus, thetotal number of solutions for the complete 2 2 nonlinearalgebraic system is s2 4.The extension of this description for the complete 3 3system and in general, for the complete n n system isstraightforward. Regarding the type of the system roots namely, real, imaginary, or complex roots - it depends ofthe values of the coefficients of the matrix A.5 ANNs as nonlinear system solversThe extension of the structure of the previously described neural network to work with the nonlinear case isstraightforward; however a few modifications are necessary.The most important of them is the usage of multiplicativeneuron types known as product units [10] in order to produce the nonlinear terms of the algebraic system (such asx2 , xy, x2 z) during training. The learning algorithm of thenetwork is the back propagation algorithm: since there isno Lyapunov energy function to be minimized, the roots of

the nonlinear system are not estimated as the outputs of theneural network, but as the weights of the synapses that jointhe input neuron with the neurons of the first hidden layer.These weights are updated during training in order to givethe desired roots. On the other hand, the weights of the remaining synapses are kept to fixed values. Some of themhave a value equal to the unity and contribute to the generation of the linear and the nonlinear terms, while someothers are set to the coefficients of the nonlinear system tobe solved. The network neurons are joined in such as way,that the total input of the output neurons - whose numbercoincides with the number of the unknowns of the system to be equal to the left hand part of the corresponding systemequation. Regarding the input layer, it has only one neuronwhose input has a always a fixed value equal to the unity.The structure of the neural network solver for the complete2 2 and 3 3 nonlinear systems are presented in the nextsections.5.1The complete 2 2 nonlinear systemThe complete nonlinear system with two equations andtwo unknowns has the general formα11 x2 α12 y 2 α13 xy α14 x α15 yα21 x2 α22 y 2 α23 xy α24 x α25 y β1β2(24)(25)and the structure of the neural network that solves it, isshown in Figure 1. From Figure 1 it is clear that the neural network is composed of eight neurons grouped in fourlayers as follows: Neuron N1 belongs to the input layer,neurons N2 and N3 belong to the first hidden layer, neuronsN4 , N5 and N6 belong to the second hidden layer, while,neurons N7 and N8 belong to the output layer. From theseneurons, the input neuron N1 gets a fixed input with a valueequal to the unity. The activation function of all the network neurons is the identity function, meaning that eachneuron copies its input to its output without modifying it this property in mathematical form is expressed as Ik Ok(k 1, 2, . . . , 8) where Ik and Ok is the total input and thetotal output of the kth neuron, respectively. Regarding thesynaptic weights they are denoted as Wij W (Ni Nj )- in other words, Wij is the weight of synapse that joins thesource neuron Ni with the target neuron Nj - and they areassigned as follows:W12 x,W13 yThese weights are updated with the back propagation algorithm and after the termination of the training operation theykeep the values (x, y) of one of the roots of the non-linearalgebraic system.ββααW25 W25 W24 W34 W36 W36 1These weights have a fixed value equal to the unity; the roleof the corresponding synapses is simply to supply to themultiplicative neurons N5 , N4 and N6 the current values ofx and y, to form the quantities x2 , xy, and y 2 , respectively.In the above notation, the superscripts α and β are used todistinguish between the two synapses that join the neuronN2 with the neuron N5 as well as the neuron N3 with theneuron N6 .W57 α11 , W67 α12 , W47 α13 , W27 α14 , W37 α15W58 α21 , W68 α22 , W48 α23 , W28 α24 , W38 α25The values of the above weights are fixed, and equal to theconstant coefficients of the nonlinear system.One of the remarkable features of the above networkstructure is the use of the product units N5 , N4 and N6 ,whose total input is not estimated as the sum of the individual inputs received from the input neurons, but as theproduct of those inputs. If we denote with Is and Ip thetotal net input of a summation and a product unit respectively, each one connected to N input neurons, then, if thisneuron is the ith neuron of some layer, the above total inputQNPNWis estimated as Ia j 1 Wij xj and Im j 1 xj ijwhere Wij is the weight of the synapse joining the jth inputneuron with the current neuron, and xj is the input comingfrom that neuron. It can be proven [10] that if the neuralnetwork is trained by using the gradient descend algorithm,then, the weight update equation for the input synapse connecting the jth input unit with the th hidden product unithas the form0hhw j(t 1) w h(t) ηflh (nethp )eρp µ¶ ln xjp cos(πξp ) πI p sin(πξp ) MX0o(dpk opk )fko (netopk )wk (26)k 1where nethp is the total net input of the product unit andfor the pth training pattern, xjp is the jth component of thatpattern, netopk is the total net input of the kth output neuron,dpk and opk are the kth component of the desired and thereal output vector respectively, M is the number of outputneurons, I p has a value of 0 or 1 according to the sign ofxjp , n is the learning rate value, while, the quantities ρp and ξp are defined asρpm NXj 1hwmjln xjp and ξpm NXhpwmjImj 1respectively. In the adopted non linear model, all the threshold values are set to zero; furthermore, since the activationfunction of all neurons is the identity function, it is clear

Figure 1. The structure of the 2 2 nonlinear system solverthat the total input calculated above is also the output thatthis neuron sends to its output processing units.The neural network model presented so far, has been designed in such a way that the total input of the output neuronN7 to be the expression α11 x2 α12 y 2 α13 xy α14 x α15 y and the total input of the output neuron N8 to be theexpression α21 x2 α22 y 2 α23 xy α24 x α25 y. To understand this, let us identify all the paths that start from theinput neuron N1 and terminate to the output neuron N7 , aswell as the input value associated with them. These pathsare the following: Path P1 : it is defined as N1 N2 N5 N7 .In this path, the neuron N2 gets a total input I2 W12 · O1 1 · x x and forwards it to the neuron N5 . There are two synapses between N2 and N5with a fixed weight equal to the unity; since N5 is amultiplicative neuron, its total input is estimated asβαO2 ) · (W25O2 ) (1 · x) · (1 · x) x2 .I5 (W25Furthermore, since neuron N5 works with the identityfunction, its output O5 I5 x2 is sent to the output neuron N7 multiplied by the weight W57 α11 .Therefore, the input to N7 emerged from the path P1is ξ17 α11 x2 . Path P2 : it is defined as N1 N3 N6 N7 .Working in a similar way, the output of N3 is O3 I3 W13 O1 y, the output of N6 is O6 I6 βα(W36O3 ) · (W36O3 ) y 2 , and the total input to N7from the path P2 is ξ27 α12 y 2 . Path P3 : it is defined as N1 N2 N4 N7 . Inthis case the total input of the multiplicative neuron N4is equal to I4 O4 W24 O2 W34 O3 xy andtherefore, the contribution of the path P3 to the totalinput of the neuron N7 is equal to ξ37 α13 xy. Path P4 : it is defined as N1 N2 N7 and contributes to N7 an input of ξ47 α14 x. Path P5 : it is defined as N1 N3 N7 and contributes to N7 an input of ξ57 α15 y.Since N7 is a simple additive neuron, its total input receivedfrom the five paths defined above is equal toI7 7Xξi7 α11 x2 α12 y 2 α13 xy α14 x α15 y (27)i 1Working in the same way, we can prove that the total inputsend to the output neuron N8 is equal toI8 8Xξi8 α21 x2 α22 y 2 α23 xy α24 x α25 y (28)i 1Therefore, if the back propagation algorithm will be usedwith the values of β1 and β2 as desired outputs, the weightsW12 and W13 will be configured during training in such away, that when the trained network gets as input the unity,the network output will be the coefficients of the secondpart of the non-linear system. But this property means thatthe values of the weights W12 and W13 will be one of theroots (x, y) of the nonlinear system - which of the roots isactually the estimated root is something that requires furtherinvestigation.

5.2The complete 3 3 nonlinear systemThe complete 3 3 nonlinear system is given by theequationsβi αi01 x3 αi02 y 3 αi03 z 3 αi04 x2 y αi05 xy 2 αi06 x2 z αi07 xz 2 αi08 y 2 z αi09 yz 2 αi10 xyz αi11 xy αi12 xz αi13 yz αi14 x2 αi15 y 2 αi16 z 2 αi17 x αi18 y αi19 z(29)(i 1, 2, 3) and the structure of the neural network used forsolving it, is shown in Figure 2. In fact, all the gray-coloredneurons representing product units belong to the same layer(their specific arrangement in Figure 2 is only for viewingpurposes) and therefore the proposed neural network structure consists of twenty three neurons grouped to four different layers with the neuron activation function to be again theidentity function. Due to the great complexity of this network it is not possible to label each synapse of the net as inFigure 1. However, the assignment of the synaptic weightsfollows the same approach. Therefore, after training, thecomponents (x, y, z) of the identified root are the weightsW12 x, W13 y and W14 z. The weights of all theinput synapses to all the product units are set to the fixedvalue of unity to generate the non linear terms of the algebraic system, while, the weights of the synapses connectingthe hidden with the three output neurons are set as W18,21W19,21W09,21W05,21W10,21W06,21 α101 α102 α103 α104 α105 α106 α107 α108 α109 α110 α111 α112 ,22W19,22W09,22W05,22W10,22W06,22 α201 α202 α203 α204 α205 α206 α207 α208 α209 α210 α211 α212 ,23W19,23W09,23W05,23W10,23W06,23 α301 α302 α303 α304 α305 α306 α307 α308 α309 α310 α311 α312 α313W07,21 α114W11,21 α115W07,22 α214W11,22 α215W07,23 α314W11,23 α315W08,21 α116W02,21 α117W08,22 α216W02,22 α217W08,23 α316W02,23 α317W03,21 α118W04,21 α119W03,22 α218W04,22 α219W03,23 α318W04,23 α319The proposed neural network shown in Figure 2 has beendesigned in such a way, that the total input Oi sent to ithneuron to be equal to the left hand part of the ith equationof the nonlinear algebraic system (i 1, 2, 3). This factmeans that if the network trained by the back propagationalgorithm and by using the constant coefficients of the nonlinear system β1 , β2 and β3 as the desired outputs, then,after training, the weights W12 , W13 and W14 will containthe values (x, y, z) of one of the roots of the nonlinear algebraic system.5.3The complete n n nonlinear systemThe common feature of the proposed neural networksolvers presented in the previous sections, is their structure,characterized by the existence of four layers: (a) an inputlayer with only one neuron whose output synaptic weightswill hold after training the components of one of the systemroots, (b) a hidden layer of summation neurons that generate the linear terms x, y in the first case and x, y, z in thesecond case, (c) a hidden layer of product units that generate the nonlinear terms by multiplying the linear quantitiesreceived from the previous layer, and (d) an output layer ofsummation neurons, that generate the left hand parts of thesystem equations and estimate the training error by comparing their outputs to the values of the associated fixed coefficients of the nonlinear system under consideration. It hasto be mentioned however that even though this structure hasbeen designed to reproduce the algebraic system currentlyused, it is not nessecarily the optimum one and the possible optimization of it, is a subject of future research. In thisproject the focus is given to the complete nonlinear systemthat contains not only the terms with a degree of s but alsoall the smaller degrees with a value of m (1 m s).The previously described structure, characterizes alsothe neural solver for the complete n n nonlinear algebraic system. Based on the fundamental theory of nonlinear algebraic systems presented in previous sections, it isnot difficult to see that for a system of n equations with nunknowns characterized by a common degree s1 s2 · · · sn s, the total number of the linear as well as ofthe nonlinear terms of any degree m (1 m s) is equaltoNn s sXm 1Mn m sX1(n m 1)!(30)(n 1)! m 1m!Table 1 presents the number of those terms for the valuesn 1, 2, 3, 4, 5 and for a degree s n, a property thatcharacterizes the nonlinear algebraic systems studied in thisresearch. From this table, one can easily verify, that for values n 2 and n 3, the associated number of the systemterms is equal to N2 2 5 and N3 3 19, in complete accordance with the number of terms of the analytic equationsfor the complete 2 2 and 3 3 systems, presented in theprevious sections.

Figure 2. The structure of the 3 3 nonlinear system solverFor any given pair (n, s), the value of Nn s identifiesthe number of the neurons that belong to the two hiddenlayers of the neural solver. In fact, the number of product units belonging to the second hidden layer is equal toP Nn s n, since, n of those neurons are the summationunits of the first hidden layer that generate the linear termsz1 , z2 , . . . , zn . Regarding the total number of neurons in thenetwork is is clearly equal to N Nn s n 1, since theinput layer contains one neuron with a constant input equalto the unity, and the output layer contains n neurons corresponding to the n fixed parameters of the nonlinear system.Therefore, the structure of the neural simulator, contains thefollowing layers: (a) an input layer with one summationunit, (b) a hidden layer of n summation units correspondingto the linear terms with a degree i 1, (c) a hidden layerwith Nn s n produc

ods for solving nonlinear systems of equations that are com-binations of the nonlinear ABS methods and quasi-Newton methods. Another interesting class of methods have been proposed by Kublanovskaya and Simonova [8] for estimat-ing the roots of m nonlinear coupled algebraic equations