Zhivko Dimitrov Georgiev
email: zhdgeorg@tu-sofia.bg

Analysis and Synthesis of Oscillator Systems Described by Perturbed Single Well Duffing Equations

Keywords: Duffing equation, Duffing oscillators, limit cycles, Melnikon function, elliptic functions, synthesis of oscillator systems

The paper presents an analysis as well as a synthesis of oscillator systems described by single well Duffing equations under polynomial perturbations of fourth degree. It is proved that such a system can have a unique hyperbolic limit cycle. There has been obtained an analytical condition for the arising of a limit cycle and an equation giving the parameters of this limit cycle. There has been proposed a method for the synthesis of oscillator systems of the considered type, having preliminarily assigned properties. The synthesis consists of an appropriate choice of the perturbation coefficients in such a way, that the oscillator equation should have a preliminary assigned limit cycle. Both the analysis and the synthesis are performed with the help of the Melnikov function.

Introduction

The paper presents an analysis as well as a synthesis of oscillator systems allowing self-sustained oscillations. The results obtained are based on the qualitative investigation of dynamic systems, whereas the self-sustained oscillations are regarded as limit cycles on the phase plane.

The Melnikov theory allows us to establish the number, position, shape and stability of the limit cycles arising in autonomous perturbed planar systems. According to this theory the problem of finding the limit cycles is replaced by an equivalent problem of finding the zeros of a given function -the Melnikov function for autonomous systems [START_REF] Andronov | Theory of Bifurcations of Dynamical Systems on the Plane[END_REF][START_REF] Bautin | Methods and Tools for Qualitative Analysis of Dynamical Systems on the Plane[END_REF][START_REF] Perko | Differential Equation and Dynamical Systems[END_REF].

Using the Melnikov theory, there has been performed an analysis of the limit cycles in oscillator systems described by single well Duffing equations under polynomial perturbations of 4 th degree. It is proved that such a system can have a unique hyperbolic (simple) limit cycle. Conditions are derived for the arising of a limit cycle (or periodic oscillations) and the establishing whether this limit cycle is stable or unstable. Moreover, a method is suggested for the synthesis of oscillator systems having preliminarily assigned properties.

In general, the synthesis of oscillators systems includes the following two stages: first, finding a differential equation to satisfy the preliminary assigned requirements and second, obtaining an electronic circuit corresponding to this differential equation. In this paper, it should be emphasized, that we are only interested in the first stage. The second stage is a separate problem and we are not examining it here. It is worth mentioning that detailed information concerning the circuit simulation of a given differential equation is presented in [START_REF] Itoh | Synthesis of Electronic Circuits for Simulating Nonlinear Dynamics[END_REF] (see also the references therein cited).In this expose, the synthesis of an oscillator system consists of choosing the perturbation coefficients in such a way, that the perturbed Duffing equation should have preliminary assigned limit cycles, respectively self-sustained oscillations.

The Melnikov theory is a successful approach for the synthesis of oscillator systems with preliminary assigned solutions and desired characteristics. Some applications of the Melnikov function for synthesizing generalized Van der Pol systems (perturbed harmonic oscillators) and some types of perturbed Duffing systems are given in References [START_REF] Savov | Using the Melnikov Function for a Synthesis of Generalized Van der Pol Systems[END_REF][START_REF] Savov | Analysis and Synthesis of Perturbed Duffing Oscillators[END_REF].

The analysis of the Melnikov function is a separate mathematical problem and there exist various mathematical techniques for its solving [START_REF] Chow | Normal Forms and Bifurcation of Planar Vector Fields[END_REF][START_REF] Christopher | Limit Cycles of Differential Equations[END_REF][START_REF] Rousseau | Bifurcation Methods in Polynomial Systems[END_REF][START_REF] Bonnin | Harmonic Balance, Melnikov Method and Nonlinear Oscillators under Resonant Perturbation[END_REF][START_REF] Sanjuan | Lienard Systems, Limit Cycles, and Melnikov Theory[END_REF][START_REF] Petrov | Non-oscillation of elliptic integrals[END_REF]. In the case of a polynomial perturbed Hamiltonian system the Melnikov function is generally an Abelian integral. Finding the zeros of an Abelian integral is known as the weakened Hilbert's 16 th problem (which still remains open [START_REF] Li | Hilbert's 16 th Problem and Bifurcations of Planar Polynomial Vector Fields[END_REF][START_REF] Han | Some Bifurcation Methods of Finding Limit Cycles[END_REF]).

In this paper expressing the unperturbed solution through the Jacobi elliptic functions allows us to compute the Melnikov function and the result is a function involving complete elliptic integrals. There exist some advantages in the study of the Melnikov function involving complete elliptic integrals compared to the study of the Melnikov function involving Abelian integrals. First of all, the complete elliptic integrals have been extensively studied -their properties and series expansions are well known. Furthermore, there exist a lot of handbooks, tables and computer software giving the values of the complete elliptic integrals [START_REF]Handbook of Mathematical Functions[END_REF][START_REF] Byrd | Handbook of Elliptic Integrals for Engineers and Physicists[END_REF][START_REF] Dwight | Tables of Integrals and other Mathematical Data[END_REF][START_REF] Oberhettinger | Anwendung der Elliptishen Funktionen in Physik und Technik[END_REF][START_REF] Janke | Tafeln Höherer Funktionen[END_REF][START_REF] Lawden | Elliptic Functions and Applications[END_REF]. All this facilitates the engineering applications of the obtained results. In addition, the analytical expression of the interdependence between both the Hamiltonian level h and the modulus of the complete elliptic integrals k is an essential element of the synthesis of oscillator systems.

In order to analyze the Melnikov function, it is necessary to prove that a given function involving complete elliptic integrals is a monotonic function. The monotonicity of this function is proved with the help of appropriate Picard-Fuchs and Riccati equations. Then this monotonic function is used to establish that the Melnikov function has a single simple zero. Equations of Picard-Fuchs and Riccati are widely used to prove that a given function defined as a ratio of two Abelian integrals is monotonic [START_REF] Chow | Normal Forms and Bifurcation of Planar Vector Fields[END_REF][START_REF] Christopher | Limit Cycles of Differential Equations[END_REF][START_REF] Dumortier | Perturbations from an Elliptic Hamiltonian of Degree Four[END_REF][START_REF] Zhang | A Study of the Limit Cycles Associated with a Generalized Codimension-3 Lienard Oscillator[END_REF][START_REF] Zhang | Degenerate Bifurcation Analysis on a Parametrically and Externally Excited Mechanical System[END_REF]. In this paper these equations are obtained for the functions involving complete elliptic integrals.

Brief information concerning the Melnikov function and the limit cycles bifurcations in perturbed planar Hamiltonian system is given below. It should be mentioned, that only the case of first order perturbations, respectively first order bifurcations of limit cycles from a "centre" and a first order Melnikov function are examined.

Let's consider the following system

) , , ( ) , , ( ε ε ε
ε y x q x H y y x p y H x + ∂ ∂ - = + ∂ ∂ = , (1)
y x q are analytic in R R × 2 , ε is a small parameter, 1 << ε , t d d ) ( ≡ ⋅
and t is time. Let's denote the solution of system (1) by ( )

) ( ), ( ) , ( t t y x ε ε ψ ϕ = .
In this case the solution of the unperturbed system (at

0 = ε ) is ( ) ) ( ), ( ) , ( 0 0 t t y x ψ ϕ = .
Further on, let's assume that the unperturbed system (at 0 = ε ) has at least one equilibrium point, which is a "centre", surrounded by closed concentric trajectories parameterized by the set of constant Hamiltonian levels. Then the equation on the phase plane of a given closed trajectory ) ( 0 h Γ , having a parameter h , is expressed by

R ⊂ ∈ = Ω Γ h h y x H h , ) , ( : ) ( 0 , ( 2
)
where

Ω
[ ] dt p q h M h T ∫ - = ) ( 0 0 0 0 0 0 0 0 ) 0 , , ( ) 0 , , ( ) ( ψ ψ ϕ ϕ ψ ϕ , (3)
is called First Order Melnikov function (briefly, Melnikov function). It plays a basic role in the study of limit cycles. The zeros of the Melnikov function determine the limit cycles in the system (1) emerging from the periodic orbits of the unperturbed Hamiltonian system. The stability of the bifurcating limit cycle is determined by the derivative of the function ) (h M . More details on this subject can be found in [START_REF] Blows | Bifurcation of Limit Cycles from Center and Separatrix Cycles of Planar Analytic Systems[END_REF][START_REF] Chicone | On Bifurcation of Limit Cycles from Centers[END_REF][START_REF] Chicone | Bifurcation of Limit Cycles from Quadratic Isochrones[END_REF].

In all computations further on the complete elliptic integrals of the first and second kind K and E are regarded as functions of 2 k ( k is the modulus), i.e.

) ( 2 k K K = and ) ( 2 k E E =
. Some properties and identities concerning the Jacobi elliptic function sd and the complete elliptic integrals K and E , that have been used in the presentation, are given in an Appendix at the end of the paper for quick reference. The equations in the Appendix are marked by the letter A.

Oscillator circuits and oscillator equations

Examples of AC equivalent circuits of oscillator systems governed by perturbed Duffing equations are shown in Fig. 1. The parallel circuit, shown in Fig. 1a, consists of a resistor with conductance G , an inductive element with inductance L , a nonlinear resistor NR of the type N , having a voltage-current characteristic

) ( 1 u f i =
and a nonlinear capacitive element having a charge-voltage characteristic

3 q q u β α + =
, where q is the electric charge associated with the element, and α , β are constants, 0 > α , 0 > β . The series circuit, shown in Fig. 1b, consists of a resistor with resistance R , a capacitive element with capacitance C , a nonlinear resistor NR of the type S having a current-voltage characteristic

) ( 1 i f u =
and a nonlinear inductive element having a flux-current characteristic

3 βψ αψ + = i
, where ψ is the magnetic flux associated with the element, and here, too, α ,

β are constants, 0 > α , 0 > β .
Let's consider the parallel oscillator circuit. The characteristic of the nonlinear resistor is transformed in the following way:

) ( ) ( ) ( 2 3 1 1 q f q q f u f i = + = = β α . ( 4
) Using Kirchhoff's current law we obtain 0 ) ( d d d 1 2 = + + + ∫ q f t q u u L Gu . ( 5
)
After differentiating the equation above, there follows

[ ] 0 1 ) ( d d d d 2 2 2 = + + + u L q f Gu t t q .
Taking into account the expression for u , we find

[ ] 0 d d ) ( d d d d 3 2 3 2 2 = + + + + + q L q L t q q f q G q G q t q β α β α . ( 6
)
If we perform the substitutions

) ( ) ( 2 3 q F q f q G q G ε β α - = + + , x q = , a L = α , b L = β , Equation (6)
takes its final form

0 d d d ) ( d d d 3 2 2 = + + - bx ax t x x x F t x ε , 0 > a , 0 > b . (7)
Equation ( 7) is known as a perturbed single well Duffing equation. The small parameter ε depends on the circuit parameters.

Starting from the series circuit (Fig. 1b) and using the Kirchhoff's voltage law, we get the same equation as in [START_REF] Chow | Normal Forms and Bifurcation of Planar Vector Fields[END_REF].

Melnikov function for perturbed single well Duffing equation

In our work we assume that the function ) (x F is represented by the polynomial

4 4 3 3 2 2 1 ) ( x a x a x a x a x F + + + = . ( 8
)
In this case Equation ( 7) can be rewritten as the following perturbed Hamiltonian system . Further on, we will investigate the system (9) from a point of view of the existence of limit cycles. As in section 1 let ( )

) ( ), ( ) , ( t t y x ε ε ψ ϕ =
be the solution of the perturbed system (9).

The unperturbed system (at

0 = ε ) has a Hamiltonian 4 2 2 4 2 2 1 ) , ( x b x a y y x H + + = (10)
and a single equilibrium point ) 0 , 0 ( which is a "centre" surrounded by a continuous oneparameter family of closed trajectories, which can be expressed as follows:

h x b x a y y x H h = + + = 4 2 2 0 4 2 2 1 ) , ( : ) ( Γ , ) , 0 ( ∞ ∈ h . ( 11
)
The solution of the unperturbed system ( )

) ( ), ( ) , ( 0 0 t t y x ψ ϕ =
in the time domain and its period are expressed as follows [START_REF] Savov | An Elliptic Functions Solution to a Class of the Duffing Equations[END_REF]

⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - - - = = k t k a k b k k a t x , 2 1 sd 2 1 1 2 ) ( 2 2 2 0 ϕ , ( 12a
) ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - - - = = k t k a k t k a k b k ak t y , 2 1 nd , 2 1 cd ) 2 1 ( 1 2 ) ( 2 2 2 2 0 ψ , ( 12b
)
a k k h T T 2 2 0 0 2 1 ) ( 4 ) ( - = = K , ( 13
)
where sd , cd and nd are Jacobi elliptic functions. The modulus k and the constant h are related in the following way:

⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + - = bh a a k 4 1 2 1 2 2 , 2 2 2 2 2 ) 2 1 ( ) 1 ( k b k k a h - - = , ( 14a
) 2 1 0 2 < < k , ∞ < < h 0 . (14b)
= ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ = ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ - = ∫ ∑ ∫ ∑ = = dt b a a dt a h M T n n n T n n n 0 0 0 3 0 0 4 1 0 0 0 4 1 0 ) ( ( ) ( ϕ ϕ ϕ ψ ϕ dt k t k a k b k k a c T n n n n ∫ ∑ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - - = = 0 0 2 2 2 7 2 , 2 1 sd 2 1 1 2
, where

dz k z k b k k a c a k h M n n n n ∫ ∑ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - - - = = K 0 2 2 2 2 3 1 2 2 ) , ( sd 2 1 1 2 2 1 4 ) ( . ( 15
)
As it is seen from the last result, the Melnikov function and the arising of limit cycles do not depend on the terms with coefficients having odd indices -3 c , 5 c , 7 c . This means that the terms with even degrees in the polynomial perturbation in System (9) exert no influence on the arising of limit cycles and further we will assume that 0

4 2 = = a a
. This fact is important for the circuit realization of the oscillator equation. It is worth mentioning that similar assertions are valid for the generalized Van der Pol oscillators [START_REF] Savov | Using the Melnikov Function for a Synthesis of Generalized Van der Pol Systems[END_REF]. Finally, System (9), which we investigate, is reduced to the following system with respect to the limit cycles:

3 3 3 1 ) ( bx ax y x a x a y x - - = + + = ε . ( 16
)
Under these conditions, after quite long but straightforward calculations, using Equations (A7)-(A9) and introducing the notation 2 k m = , we obtain the final expression for the Melnikov function ( )

) , ; ( ) 2 1 ( 2 1 8 ) ( ) ( ) ( 3 1 3 a a m m m a m m h M h M B M - - = = = λ , ( 17
)
where b a = λ and ) ( 5
2 ) ( ) , ; ( 2 3 1 1 3 1 m I a m I a a a m λ + = B , (18)
[ ]

E K ) 2 1 ( ) 1 ( ) 2 1 ( 3
1 ) ( 1 m m m m I - - - - = , ( 19
) K E ) 1 )( 2 ( ) 1 ( 2 ) ( 2 2 m m m m m I - - - + - = . ( 20
)
Let's remind that here and in the equations below, the complete elliptic integrals K and E are regarded as functions of m , where

2 k m =
and k is a modulus, i.e. ) (m K K = and ) (m E E = [START_REF]Handbook of Mathematical Functions[END_REF]. According to the Melnikov theory the zeros of ) (h M determine the limit cycles emerging from the periodic orbits of the unperturbed Hamiltonian systems. The stability of the limit cycle is determined by the sign of the quantity

) ( dh dM ε . Since ) (m h h = and ) (h m m = are single-valued functions the zeros of ) (h M and ) (m M coincide. In other words, if 0 m is a zero of ) (m M , i.e. 0 ) ( 0 = m M , then ) ( 0 0 m h h = is a zero of ) (h M , i.e. 0 ) ( 0 = h M . Moreover, from the relations h m m m h m h M d d d ) ( d d ) ( d d d M M = = , 0 d d > h m ,
m I a m I a m m B λ + = = B . ( 21
)
Therefore we will look for the zeros of the function

) (m B
, or the roots of the equation 0

) ( = m B
, and the signs of the derivative ) (m B′ evaluated at these roots. For further analysis it is convenient to rewrite the function ) (m B in the following way

⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + = 3 2 5 2 ) ( ) ( ) ( a m Y m I m B λ , (22)
where

[ ] K E E K ) 1 )( 2 ( ) 1 ( 2 ) 2 1 ( ) 1 ( ) 2 1 ( 3 1 ) ( ) ( ) ( 2 2 1 m m m m m m m m I m I m Y - - - + - - - - - = = . ( 23
)
The analysis of the Melnikov function (respectively the analysis of function

) (m B
), requires some preliminary results that are given in the following section.

Functions involving complete elliptic integrals and equations of Picard-Fuchs and Riccati

In order to find the zeros of the Melnikov function we need to analyze the functions . The properties of these functions, which will be needed later, are summed up in the following lemma:

Lemma 1: The following assertions hold: (a) 0 ) 0 ( 1 = I , 0 ) 2 1 ( 1 = I , 4 ) 0 ( 1 π = ′ I and the function ) ( 1 m I is positive on ) 2 1 , 0 ( , i.e. 0 ) ( 1 > m I for ) 2 1 , 0 ( ∈ m . ( 24
) (b) 0 ) 0 ( 2 = I , 0 ) 2 1 ( 2 > I , 0 ) 0 ( 2 = ′ I and
> m I , 0 ) ( 2 > ′ m I for ) 2 1 , 0 ( ∈ m . ( 25
) (c) The functions ) ( 1 m I and ) ( 2 m I satisfy the following system ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ - - + - = ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ - - 2 1 2 2 2 1 0 ) 1 ( 45 ) 2 1 ( ) 6 6 1 ( 6 d d ) 2 1 )( 1 ( 6 I I m m m m m I I m m m m , ( 26
)
which is called the Picard-Fuchs equation.

Proof: (a), (b). In the first place, using Equations (A12),we calculate the derivatives

[ ] K E ) 10 7 ( ) 2 1 ( 10 6 1 ) ( 1 m m m I - - - = ′ , (27)
[ ]

) 2 1 ( ) ( 2
15 ) 2 1 ( ) 1 ( 2 5 ) ( 1 2 m m I m m m I - = - - - = ′ E K . ( 28
)
Then the values ) 0 (

1 I , ) 2 1 ( 1 I , ) 0 ( 1 I ′ , ) 0 ( 2 I , ) 0 ( 2 I ′
are obtained by direct calculations. Moreover, using elementary properties of the complete elliptic integrals, we find

[ ] [ ] 0 ) )( 1 ( ) 2 1 ( 3
1 ) 2 1 ( ) 1 ( ) 2 1 ( 3
1 ) ( 1 > + - - - = - - - - = E E K E K m m m m m m m I , From here it follows that 0 ) ( 2 > ′ m I , i.e. the function ) ( 2 m I is monotone increasing. Further, in view of Equation (A11) we get 0 ) 2 1 ( 2 1 ) 2 1 ( 2 3 ) 2 1 ( 2 > ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ - = K E I .
From the relations 0 ) 0 ( is not a rational function, we make equal the expressions in front of K , respectively E , on the left and right sides of the system equations. In this way we obtain the coefficients A , B , C , D and this yields the system [START_REF] Chicone | Bifurcation of Limit Cycles from Quadratic Isochrones[END_REF].■

The following lemma gives the properties of the function ) (m Y .

Lemma 2:

The following assertion is valid:

(a) The function ) (m Y satisfies the following Riccati equation 2 2 2 ) 2 1 ( ) 6 6 1 ( 6 ) 1 ( 45 ) 2 1 )( 1 ( 6 m Y m m Y m m Y m m m - + + - + - - = ′ - - . ( 29
) (b) The function ) (m Y satisfies the following second order differential equation Y Y m m m Y m Y m Y m m m ′ - - - - - - - - = ′ ′ - - ) 1 ( 90 ) 2 1 ( 4 ) 2 1 ( 36 ) 2 1 ( 45 ) 2 1 )( 1 ( 6 2 . (30) (c) +∞ = = + → ) ( lim ) 0 ( 0 m Y Y m , 0 ) 2 1 ( = Y and 0 ) ( > m Y for ) 2 1 , 0 ( ∈ m . ( 31
) (d) The function ) (m Y is strictly monotone decreasing on the interval ) 2 1 , 0 ( , i.e. 0 ) ( < ′ m Y for ) 2 1 , 0 ( ∈ m . ( 32
)
Proof: (a) Differentiating ( 23) and introducing

) 2 1 )( 1 ( 6 m m m c - - = we obtain 2 2 2 1 2 2 2 1 2 1 I I c Y I I c I I I I I c Y c ′ - ′ = ′ - ′ = ′ .
Replacing in this equation 1 I c ′ and 2 I c ′ from ( 26), and making some transformations we get (29). (b) Equation ( 30) is obtained after differentiation of (29).

(c) The values of ) (m Y with + → 0 m , and ) 2 1
( Y are obtained by using l'Hôpital rule and direct calculations using the results of Lemma 1. We obtain

+∞ = ′ ′ = = + → + → ) ( ) ( lim ) ( ) ( lim ) 0 ( 2 1 0 2 1 0 m I m I m I m I Y m m , [ ]
0 ) 2 1 ( ) 2 1 ( ) 2 1 ( ) 2 3 ( 0 ) 2 1 ( ) 2 1 ( ) 2 1 ( 2 1 = - = = K E I I Y .
Moreover

) ( = ′ m Y , 0 ) ( > ′ ′ m Y for ) 2 1 , 0 ( ∈ m .
Taking into account the last relations, Equation (30) yields

0 ) 2 1 ( 4 ) ( ) 2 1 ( 36 ) ( ) 2 1 ( 45 ) ( ) 2 1 )( 1 ( 6 2 < - - - - - - = ′ ′ - - m m Y m m Y m m Y m m m .
The last result gives the inequality 0

) ( < ′ ′ m Y
, which means that the stationary point m is a point of maximum (it could be stated that all stationary points are points of maximum). This fact is in contradiction to our supposition and this proves part (d) of the lemma.■

The graph of the function ) (m Y is shown in Figure 2.

Analysis of the Melnikov Function and limit cycles

We are now able to formulate the following theorem which gives the main results concerning the limit cycles analysis in perturbed single well Duffing oscillators.

Theorem:

The system (16) (with

1 1 = a
) can have a unique hyperbolic limit cycle. The condition for the arising of a limit cycle is the fulfillment of the following inequality

0 3 < a . ( 33
)
The parameters of the limit cycle are determined by the equation

0 5 2 ) ( 3 = + a m Y λ . ( 34
)
In addition, the limit cycle is stable when 0 > ε and unstable when 0 < ε .

Proof:

As it has been emphasized above, the zeros of the Melnikov function are determined by the solution of the following equation

0 5 2 ) ( ) ( ) ( 3 2 = ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + = a m Y m I m B λ . ( 35
)
We

) ( 0 = m B , 0 ) ( ) ( ) ( 0 0 2 0 < ′ = ′ m Y m I m B . (36)
This, according to the Melnikov theory, means that for a sufficiently small 0 ≠ ε , there exists such an ε h in an

) (ε O neighbourhood of ) ( 0 0 m h h =
(see Equations (14a)), that the perturbed Hamiltonian system (16) has a hyperbolic (simple) limit cycle

) ( ε ε Γ h . The limit cycle ) ( ε ε h Γ is localized in a ) (ε O neighbourhood of the curve ) ( 0 0 h Γ and tends to ) ( 0 0 h Γ with 0 → ε . The

Synthesis of perturbed single well Duffing oscillators

The presented theory can be used to perform a synthesis of the considered oscillators. The relations between the Hamiltonian level h and the elliptic modulus k and 2 k m = given in Equation (14a) allow us to compute the value of 0 m corresponding to the preliminarily assigned Hamiltonian level 0 h . Then Equation (34), for a given 0 m , yields the coefficient 3 a in such a way that the system ( 16) is to have a limit cycle localized in a ) (ε O neighbourhood of the assigned Hamiltonian level. The synthesis procedure includes the following steps:

1) Formulation of the synthesis problem: Find an oscillator system of the type ( 16) (with

1 1 = a ), having a limit cycle localized in an ) (ε O neighborhood of the curve 0 4 2 2 0 0 4 2 2 1 : ) ( h x b x a y h Γ = + + , (37)
where the values of a , b , 0 h are given. We note that the preliminarily assigned values of a , b ( b a = λ ) and 0 h form the desired characteristics of the oscillator. In other words, these values determine the shape, position and period of the limit cycle which will arise in system [START_REF] Byrd | Handbook of Elliptic Integrals for Engineers and Physicists[END_REF]. The perturbation coefficient 3 a obtained in this way ensures that the system ( 16) is to have a simple limit cycle ) ( ε ε h Γ , which is localized in a ) (ε O neighborhood of the curve (37). The limit cycle ) ( ε ε h Γ is stable with 0 > ε , and unstable with 0 < ε . Finally, as an application of the proposed method, we shall briefly consider an example of a synthesis of an oscillator system. The computation of the numerical data is given after the synthesis procedure.

[ ]

k m = are ) 1 ( 2 ) ( ) 1 ( ) ( d )

  account Equation[START_REF] Li | Hilbert's 16 th Problem and Bifurcations of Planar Polynomial Vector Fields[END_REF] and the properties of the function sd given in the Appendix (Equations (A1) -(A3)), we see that the integrals containing

  In order to obtain the system[START_REF] Chicone | Bifurcation of Limit Cycles from Quadratic Isochrones[END_REF], we use the method of undetermined coefficients. Let's consider the system A , B , C and D are unknown expressions. Replacing the expressions 1 I , 2 I , 1 I ′ , 2 I ′ and using the fact that E K

2 ) 4 )

 24 Computation of the values of 0 m and modulus 0 k of the complete elliptic integrals: Computation of the perturbation coefficient 3 a :

Figure 1 .

 1 Figure 1. Oscillator circuits governed by perturbed Duffing equations Figure 2. Graphical plot of the function ) (m Y Figure 3. Phase portrait of system (41)

  Figure 2.

Figure 3 .

 3 Figure 3.

  is an open set of values h , for which the trajectories are not degenerated in a point

	or a separatrix loop. Let's denote the period of the unperturbed solution	(	x	,	y	)	=	( ϕ	0	( t	), ψ	0	( t	)	)	,
	corresponding to the closed trajectory	Γ	( 0 h	)	, by	( 0 h T	)	. Using these notions the following
	function

  stability of the limit cycle is determined by the sign of the quantity

										( 0 m B′ ε	)
	and it follows from Equation (36) that the limit cycle	Γ	( ε ε h	)	is stable when	ε	>	0	, and
	unstable when	ε	<	0	. The theorem is proved.■
	The limit cycle	Γ	( ε ε h	)	in practice coincides with the curve	Γ	( 0 0 h	)	and its first order
	approximation in the time domain is expressed by the functions (12).

  Equation (A7) is given in[START_REF] Byrd | Handbook of Elliptic Integrals for Engineers and Physicists[END_REF]. Equations (A8) and (A9) are obtained with the help of the recurrence formula (A4). ■

	K 0 ∫	sd	6	d	z	=	15 k	6	1 ( 1 -	k	2	)	3	(	23 k	4	-	23 k	2	+	) 8	E	+	15 ( k	6	-	34	k	4	+	27	k	2	-	) 8	K	.	(A9)
	Proof: The following relation
				E	(	m	)	K	1 (	-	m	)	+	E	1 (	-	m	)	K	(	m	)		-	K	(	m	)	K	1 (	-	m	)	=	π	2	.	(A10)
	is known as the Legendre's Relation [15-16]. At	m	= k	2 =	1	2	Equation (A10) becomes
																K	1 (	2	[ 2 )	E	1 (	2	)	-K	1 (	2	)	] 2 π =	,
	from where it follows
															2	E	1 (		2	)	-K	1 (	2	)	>	0	.	(A11)
	The derivatives of the functions	K	(m	)	and	E	(m	)	with respect to	2

Example: Find an oscillator system of the type [START_REF] Byrd | Handbook of Elliptic Integrals for Engineers and Physicists[END_REF]

. In this way we find the following system (41

The system (41) has a hyperbolic limit cycle

Γ

. The self-sustained oscillations generated by system (41) are expressed practically in the time domain as follows

The phase portrait of the system (41) for 25 . 0 = ε and initial data ) 5 . 0 , 0 ( and ) 4 , 0 ( obtained by a numerical integration is shown in Fig. 3. The numerical perfectly computations confirm the analytical results.

Conclusions

An analysis of oscillator systems described by single well Duffing equation under polynomial perturbations of 4 th degree has been presented. It has been proved that the considered system can have a unique hyperbolic limit cycle. The obtained inequality (33) represents a condition for the arising of a limit cycle and the equation (34) determines the parameters of this limit cycle.

Moreover, a method for the synthesis of oscillator systems of the considered type is proposed. The synthesis consists of finding such appropriate perturbation coefficients that the prescribed properties are fulfilled. Both the analysis and the synthesis are performed with the help of the Melnikov function.

An example of the synthesis of oscillator system is given. Together with the numerical simulation, it illustrates and confirms the presented theory.

The results obtained can be used for more complete analysis and synthesis of oscillator systems.

Appendix

In this section we give some properties and identities concerning the Jacobi elliptic function ) , sd( sd sd k z z = = and the complete elliptic integrals of the first and second kind, K and E , which are used in the presentation. Recall that the complete elliptic integrals K and E are regarded as functions of 2 k , where k is a modulus, i.e.

is a periodic and odd function with a period K 4 . Moreover the following identities are valid

The results concerning the integration of Jacobi elliptic function sd are given in two lemmas.

Lemma A1:

The following recurrence formula is valid

where n is an even number and 4 ≥ n , (

).

Proof: In order to prove the recurrence formula (A4) we start from the identity

where cd and nd are Jacobi elliptic functions. This identity can be easily obtained using the usual dependencies between the Jacobi elliptic functions [START_REF]Handbook of Mathematical Functions[END_REF][START_REF] Byrd | Handbook of Elliptic Integrals for Engineers and Physicists[END_REF][START_REF] Dwight | Tables of Integrals and other Mathematical Data[END_REF][START_REF] Oberhettinger | Anwendung der Elliptishen Funktionen in Physik und Technik[END_REF][START_REF] Janke | Tafeln Höherer Funktionen[END_REF][START_REF] Lawden | Elliptic Functions and Applications[END_REF]. Integrating (A5) and taking into account the relations (A1), we obtain


Cem Güneri
email: guneri@sabanciuniv.edu

Ferruh Özbudak
email: ozbudak@metu.edu.tr

A Relation Between Quasi-Cyclic Codes and 2-D Cyclic Codes ⋆

We consider a q-ary quasi-cyclic code C of length mℓ and index ℓ, where both m and ℓ are relatively prime to q. If the constituents of C are cyclic codes, we show that C can also be viewed as a 2-D cyclic code of size m × ℓ over Fq. If we further assume that m and ℓ are also coprime to each other, then we easily observe that the code C must be equivalent to a cyclic code. The last fact was proved earlier by Lim using a different approach.

Introduction and Background

Both 2-D cyclic codes and quasi-cyclic (q.-c.) codes are generalizations of cyclic codes. We introduce definitions and basic facts about these codes in this section. Although there are numerous useful references on both classes of codes, we refer to [START_REF] Güneri | Artin-Schreier curves and weights of two-dimensional cyclic codes[END_REF][START_REF] Ling | On the algebraic structure of quasi-cyclic codes I: finite fields[END_REF] since our presentation will follow the notation and approach in these articles.

Let F q be the finite field with q elements. We assume throughout that m, ℓ > 1 are integers both of which are relatively prime to q. Let C be a linear code over F q of length mℓ whose codewords are viewed as m × ℓ arrays, i.e. c ∈ C is written as

c =     
c 00 , . . . , c 0,ℓ-1 c 10 , . . . , c 1,ℓ-1 . . .

c m-1,0 , . . . , c m-1,ℓ-1      . ( 1
)
If C is closed under row and column shift of codewords, then we call C a 2-D cyclic code of size m × ℓ. As in the case of cyclic codes, we have a polynomial ⋆ The authors were partially supported by T ÜB İTAK under Grant No. TBAG-109T672.

(ideal) representation for 2-D cyclic codes. For this, consider the following F qspace isomorphism

F m×ℓ q -→ F q [x, y]/(x m -1, y ℓ -1) a ij -→ m-1 i=0 ℓ-1 j=0 a ij x i y j , (2)
where F m×ℓ q denotes the space of all m × ℓ arrays. Via this identification we can conclude that a linear subspace C of F m×ℓ q is 2-D cyclic if and only if C is an ideal when viewed as a subset of F q [x, y]/(x m -1, y ℓ -1). Note that we can also think of a 2-D cyclic code as an ideal of the group algebra

F q [Z m × Z ℓ ].
Let η be a primitive ℓ th root of unity and ξ be a primitive m th of root unity in the rest of this document. We take both of these elements in some extension F of F q . Consider the following set:

Ω = {(ξ i , η j ); 0 ≤ i ≤ m -1, 0 ≤ j ≤ ℓ -1}. ( 3
)
We define the F q -conjugacy class containing (ξ i , η j ) to be

[(ξ i , η j )] = (ξ i , η j ), (ξ iq , η jq ), . . . , (ξ iq δ-1 , η jq δ-1 ) , (4)
where δ is the least common multiple of the degrees of ξ i and η j over F q . Since a 2-D cyclic code C is an ideal of F q [x, y]/(x m -1, y ℓ -1), we can write it as J/(x m -1, y ℓ -1) for some ideal J of F q [x, y]. We define the zero set Z(C) of C as the set of zeros of the ideal J over the algebraic closure Fq . Note that an F q -conjugacy class of an element (ξ i , η j ) ∈ Ω can be described by any member of the class [(ξ i , η j )]. The set formed by picking exactly one representative from each conjugacy class in Z(C) is called a basic zero set for C and denoted by BZ(C). Since J contains the ideal (x m -1, y ℓ -1), its zeros lie in the set Ω. It can be shown that Z(C) is a disjoint union of some F q -conjugacy classes of elements in Ω. Moreover, Z(C) (hence BZ(C)) uniquely determines the code C and dim(C) = |Ω \ Z(C)|. It is easy to observe that the dual C ⊥ (with respect to the usual Euclidean inner product) of a 2-D cyclic code is also 2-D cyclic. The zero set of C ⊥ is related to that of C by the formula

Z(C ⊥ ) = Ω \ Z(C) -1 , where Z(C) -1 := (ξ -i , η -j ) : (ξ i , η j ) ∈ Z(C) .
A linear code C ⊂ F mℓ q is called quasi-cyclic of index ℓ if it is invariant under shift of codewords by ℓ units. Note that a q.-c. code of index 1 is nothing but a usual cyclic code. We note that our notation will be the same as the notation of [START_REF] Ling | On the algebraic structure of quasi-cyclic codes I: finite fields[END_REF].

Remark 1. It is advantageous to think of length mℓ codewords of a q.-c. as m × ℓ arrays as in (1) since one immediately realizes that being closed under shift by ℓ units amounts to being closed under row shift. Hence, we also see clearly that any 2-D cyclic code of size m × ℓ is a q.-c. code of length mℓ and index ℓ.

Define the ring R := F q [Y ]/(Y m -1) and recall that a linear code of length ℓ over R is nothing but an R-submodule of R ℓ . Consider the map φ :

F mℓ q -→ R ℓ c =      c 00 , . . . , c 0,ℓ-1 c 10 , . . . , c 1,ℓ-1 . . . c m-1,0 , . . . , c m-1,ℓ-1      -→ (c 0 (Y ), c 1 (Y ), . . . , c ℓ-1 (Y )), (5)
where

c j (Y ) := m-1 i=0 c ij Y i = c 0j + c 1j Y + c 2j Y 2 + • • • + c m-1,j Y m-1 ∈ R for each 0 ≤ j ≤ ℓ -1.
In other words, we construct an element of R from each column of c. The map φ induces a one-to-one correspondence between index ℓ q.-c. codes of length mℓ over F q and linear codes of length ℓ over R. For the case ℓ = 1, this amounts to the well-known polynomial (ideal) presentation of cyclic codes in R.

Next, we try to understand the structure of the ring R better. Since m is relatively prime to q, polynomial Y m -1 is separable. Let

Y m -1 = f 1 f 2 . . . f r (6)
be its unique factorization into irreducible polynomials in F q [Y ]. Taking reciprocal polynomials of both sides above we obtain

Y m -1 = -f * 1 f * 2 . . . f * r , (7)
where f * i denotes the reciprocal of f i . Using the fact that the reciprocal of an irreducible polynomial is also irreducible, we see that equations ( 6) and (7) are two factorizations of the same polynomial into irreducibles. We rewrite Y m -1 as

Y m -1 = δg 1 g 2 . . . g s h 1 h * 1 . . . h t h * t , (8)
where g i 's are those f j 's which are associate to their own reciprocals and h i , h * i 's are the remaining f j 's grouped in pairs (so, 2t + s = r). Note that for a primitive m th root of unity ξ, the roots of irreducible factors in (8) are powers of ξ.

By Chinese Remainder Theorem we have

R = s i=1 F q [Y ]/(g i ) ⊕   t j=1 F q [Y ]/(h j ) ⊕ F q [Y ]/(h * j )   . ( 9
)
We will use the following notation for the direct factors above:

G i = F q [Y ]/(g i ) H ′ j = F q [Y ]/(h j ) H ′′ j = F q [Y ]/(h * j )
Since the polynomials are irreducible, each of the quotients above are finite fields (e.g. [G i : F q ] = deg g i ).

Note that (9) implies

R ℓ = s i=1 G ℓ i ⊕   t j=1 (H ′ j ) ℓ ⊕ (H ′′ j ) ℓ   . (10)
Hence, any linear code C of length ℓ over R can be decomposed as

C = s i=1 C i ⊕   t j=1 C ′ j ⊕ C ′′ j   , (11)
where

C i , C ′ j , C
′′ j are linear codes of length ℓ over the fields G i , H ′ j , H ′′ j , respectively. These length ℓ linear codes over various extensions of F q are called the constituents of C.

In [2, Theorem 4.3], a trace representation for 2-D cyclic codes is obtained starting with a basic zero set of the dual code. In [5, Theorem 5.1], a trace representation for q.-c. codes is obtained based on the decomposition in (9) and the constituent codes. Both results are extensions of the following well-known trace representation of cyclic codes.

Theorem 1. ([6, Proposition 2.1])

Let n be relatively prime to q and β be a primitive n th root of unity in some extension E of F q . Let C be a q-ary cyclic code of length n whose dual's basic zero set is BZ(C ⊥ ) = β i1 , β i2 , . . . , β is (i.e. the generator polynomial of C ⊥ is j m j (x) where m j (x) ∈ F q [x] is the minimal polynomial of β ij over F q ). Then we have

C = Tr E/Fq (c 1 β ti1 + • • • + c s β tis ) 0≤t≤n-1 : c 1 , . . . , c s ∈ E . ( 12
)
2 Results

We use the trace representations in [2, Theorem 4.3] and [START_REF] Ling | On the algebraic structure of quasi-cyclic codes I: finite fields[END_REF]Theorem 5.1] to obtain the two results presented in this section. In the following, we investigate a q.-c. code whose nonzero constituents are full ambient spaces.

Theorem 2. Suppose C is a q.-c. code over F q of length mℓ and index ℓ of the form

C = G ℓ i1 ⊕ • • • ⊕ G ℓ ie (H ′ j1 ) ℓ ⊕ • • • ⊕ (H ′ j f ) ℓ (H ′′ k1 ) ℓ ⊕ • • • ⊕ (H ′′ k h ) ℓ , (13
) where {i 1 , . . . , i e } ⊂ {1, . . . , s} and {j 1 , . . . , j f }, {k 1 , . . . , k h } ⊂ {1, . . . , t} (cf. Eqns (10) and ( 11)). Let u i1 , . . . , u ie , v j1 , . . . , v j f , w k1 , . . . , w k h denote fixed representatives of q-cyclotomic cosets mod m corresponding to the fields appearing in (13) (i.e. corresponding to the powers of ξ among the roots of g i , h j , h * j , respectively). Then

dim(C) = ℓ dim(D) and d(C) = ℓd(D),
where D is the q-ary cyclic code of length m whose dual's basic zero set is BZ(D ⊥ ) = {ξ -ui 1 , . . . , ξ -ui e , ξ -vj 1 , . . . , ξ -vj f , ξ -w k 1 , . . . , ξ -w k h }.

Recall from Remark 1 that any 2-D cyclic code of size m × ℓ is a q.-c. code of length mℓ and index ℓ. The next result shows that the converse of this statement holds for q.-c. codes with cyclic constituents in (11). Namely, let C be q.-c. of length mℓ, index ℓ over F q such that

C = (C 1 ⊕ • • • ⊕ C e ) (D 1 ⊕ • • • ⊕ D f ) (E 1 ⊕ • • • ⊕ E h ) , (14)
where

C a ⊂ G ℓ ia (1 ≤ a ≤ e), D b ⊂ (H ′ j b ) ℓ (1 ≤ b ≤ f ), E d ⊂ (H ′′ k d ) ℓ (1 ≤ d ≤ h
) are all cyclic codes of length ℓ over various extensions of F q . We define these cyclic constituent codes by their duals' basic zero sets:

BZ(C ⊥ 1 ) = n η x 1,1 , . . . , η x 1,k 1 o BZ(D ⊥ 1 ) = n η y 1,1 , . . . , η y 1,m 1 o BZ(E ⊥ 1 ) = n η z 1,1 , . . . , η z 1,n 1 o . . . . . . . . . BZ(C ⊥ e ) = n η x e,1 , . . . , η x e,ke o BZ(D ⊥ f ) =  η y f,1 , . . . , η y f,m f ff BZ(E ⊥ h ) = n η z h,1 , . . . , η z h,n h o (15)
As in Theorem 2, let u i1 , . . . , u ie , v j1 , . . . , v j f , w k1 , . . . , w k h denote fixed representatives of q-cyclotomic cosets mod m corresponding to the fields G ia , H ′ j b , H ′′ k d , respectively. Then we have: Theorem 3. Let C be a q.-c. code over F q of length mℓ and index ℓ such that

C = (C 1 ⊕ • • • ⊕ C e ) (D 1 ⊕ • • • ⊕ D f ) (E 1 ⊕ • • • ⊕ E h ) ,
where each constituent is a cyclic code of length ℓ over various extensions of F q . Assume that the cyclic constituents are defined by their duals' basic zero sets as in (15). Then C is a 2-D cyclic code of size m × ℓ over F q whose dual's basic zero set is

BZ(C ⊥ ) = S 1 ∪ • • • ∪ S e ∪ S ′ 1 ∪ • • • ∪ S ′ f ∪ S ′′ 1 ∪ • • • ∪ S ′′ h , where S1 = ˘(ξ -u 1 , η x 1,1 ), . . . , (ξ -u 1 , η x 1,k 1 ) ¯S′ 1 = ˘(ξ -v 1 , η y 1,1 ), . . . , (ξ -v 1 , η y 1,m 1 ) . . . . . . Se = ˘(ξ -ue , η x e,1 ), . . . , (ξ -ue , η x e,ke ) ¯S′ f = ˘(ξ -v f , η y f,1 ), . . . , (ξ -v f , η y f,m f ) S′′ 1 = ˘(ξ -w 1 , η z 1,1 ), . . . , (ξ -w 1 , η z 1,n 1 ) . . . S ′′ h = ˘(ξ -w h , η z h,1 ), . . . , (ξ -w h , η z h,n h ) ¯(16)
Remark 2. Note that we have taken both m, ℓ to be relatively prime to q so far. Let us make a further assumption that m and ℓ are also relatively prime to each other. Recall that we can view a 2-D cyclic code C of size m × ℓ as an ideal in the group algebra F q [Z m × Z ℓ ]. With the new assumption, we have an isomorphism of groups Z m × Z ℓ ∼ = Z mℓ . Hence, C can be viewed as an ideal in the group algebra F q [Z mℓ ] ∼ = F q [T ]/(T mℓ -1). This ring isomorphism just changes the positions of coefficients since it is induced from an isomorphism of the groups involved. Hence, if we have a q.-c. code C of length mℓ, index ℓ such that m, ℓ are both relatively prime to q and also coprime to each other, and if the constituents of C are cyclic codes, then C is equivalent to a usual cyclic code of length mℓ. This was observed by Lim in [START_REF] Lim | Quasi-cyclic codes with cyclic constituent codes[END_REF]. A special case of this remark (for a binary q.-c. code with cyclic constituents of length 5ℓ, index ℓ with ℓ relatively prime to both 5 and 2) was also observed by Bracco et al in [START_REF] Bracco | On quintic quasi-cyclic codes[END_REF].

Remark 3. General weight estimates for 2-D cyclic codes have been obtained in [START_REF] Güneri | Artin-Schreier curves and weights of two-dimensional cyclic codes[END_REF] and [START_REF] Güneri | Multidimensional cyclic codes and Artin-Schreier type hypersurfaces over finite fields[END_REF]. Realizing some q.-c. codes as 2-D cyclic codes makes it possible to use these estimates in the context of q.-c. codes.


Sergey Bezzateev

Natalia Shekhunova

Class of Binary Generalized Goppa Codes Perfect in Weighted Hamming Metric

The class of the binary generalized Goppa codes is offered. It is shown that the codes of this class are on the Hamming bound constructed for a weighted metric.

Introduction

It is well-known that the problem of construction of good code for a given channel consists of two parts: a choice of a metric agreed with channel errors and search of a set of vectors with given metric properties as a correcting code. A great number of results of algebraic coding theory are codes for channels whose errors are agreed with the Hamming metric (HM). However, there are channels for which the error probability depends on a codeword position [START_REF] Litwin | The Principles of OFDM[END_REF][START_REF] Chou | Interleaving and error correction in volume holographic memory systems[END_REF][START_REF] Sehlstedt | Nonuniform Bit-Energy Distribution in Wireless Video Frames[END_REF][START_REF] Hossein | Nonuniform Error Correcting Using Low-Density Parity-Check Codes[END_REF]. It means that the distribution of errors along codeword length is nonuniform. These channels are described by the model with nonuniform distribution of errors. Moreover, there are channels in which it is necessary to correct all errors from a certain set containing vectors with various configuration and weight [START_REF] Park | Error-Pattern-Correcting Cyclic Codes Tailored to a Prescribed Set of Error Cluster Patterns[END_REF][START_REF] Deza | Codes Correcting an Arbitrary Sets of Errors[END_REF]. The weighted Hamming metric (WHM) is a metric agreed with errors of such channels [START_REF] Deza | Encyclopedia of Distances[END_REF].

The paper presents a class of binary generalized Goppa codes constructed for the weighted Hamming metric (WHM). It is shown that parameters of the proposed class of codes are on the Hamming bound constructed for the WHM.

2 Weighted Hamming Metric Definition 1. [START_REF] Deza | Encyclopedia of Distances[END_REF] The distance between vectors a = (a 1 a 2 . . . a n ) and b = (b 1 b 2 . . . b n ) in the WHM is defined by a function d wH (a, b) that is given as

d wH (a, b) = n i=1 w i • d(a i , b i ) , where w i > 0, d(a i , b i ) = 1 , if a i = b i and d(a i , b i ) = 0 , if a i = b i .
The value w i and vector w = (w 1 , w 2 , . . . , w n ) are called as a position i weight and a vector of weights of positions, respectively. Definition 2. A set v = {v 1 , v 2 , . . . , v l } , l ≤ n consisting of all the different values of vector w components is called as a set of weights of the weighted Hamming metric.

The length of the vector a can be written as n = n 1 + n 2 + . . . + n l , where n i is a number of vector a positions with a weight v i .

Let us consider a sphere of radius τ with a center in a vector a = (a 1 a 2 . . . a n ) in a n-dimensional binary vector space V with the WHM given in it. It is clear that this sphere contains all the vectors b = (b 1 b 2 . . . b n ), such that

d wH (a, b) ≤ τ.
Then the volume of the sphere W τ n is determined by

W τ n = {τ1,...,τ l }: l P i=1 τivi≤τ l i=1 n i τ i .
The Hamming bound in the WHM can be written in the following way:

M • W τ n ≤ 2 n ,
where n = n 1 +n 2 +. . .+n l is a length of codeword, M is a number of codewords, d wH is a minimum distance of the code in the WHM and τ = ⌊

M • W τ n = 2 n .
3 Special class of binary generalized Γ (L, G) codes Generalized Γ (L, G) codes [START_REF] Shekhunova | Cyclic (L, G)-codes[END_REF] are natural generalization of the classical Γ (L, G) codes introduced by V. Goppa [START_REF] Goppa | A new class of linear error correcting codes[END_REF]. They are determined by a Goppa polynomial G(x) and set of numerators of codeword positions L.

In classical Γ (L, G) codes, rational fraction functions

f i (x) = 1 x-αi , G(α i ) = 0 are elements of the set L.
In generalized Γ (L, G) codes, rational fraction functions f i (x) = νi(x) ui(x) , where

u i (x) is a polynomial over GF (2 m ) , deg u i (x) > 1 and (u i (x), G(x)) = 1, (u i (x), u j (x)) = 1, deg ν i (x) < deg u i (x) < deg G(x), (u i (x), ν i (x))) = 1
for all i, j = 1, . . . , n , i = j are elements of the set L. It should be noted that for binary Γ (L, G) codes numerators of functions f i (x) are formal derivatives from denominators for both classical and generalized codes.

In this paper, we consider binary generalized Γ (L, G) codes with an irreducible Goppa polynomial G(x) and set L =

u ′ i (x) ui(x)
, where u i (x) is a normalized irreducible polynomial over the field GF (2 m ), deg u i (x) ≤ l and u ′ i (x) is a formal derivative of polynomial u i (x) . Let us use the following Lemma in order to set the special class of codes that is considered here. Lemma 1. For any separable polynomial ϕ(x) = i u i (x) with coefficients from GF (2 m ) the following relation holds:

ϕ ′ (x) ϕ(x) = i u ′ i (x) u i (x)
,

where ϕ ′ (x), u ′ i (x) are formal derivatives of ϕ(x) and u i (x) respectively.

Now let us give a definition of a generalized binary Γ (L, G) -code agreed with the WHM .

Definition 4. [START_REF] Bezzateev | Generalized Goppa codes for correcting localized errors[END_REF][START_REF] Bezzateev | One generalization of Goppa codes[END_REF] The set of all binary vectors a = (a

(1) 2 . . . a (1) n1 a

(2)

1 . . . a (l) 1 . . . a (l) n l ) of length n = n 1 + n 2 + . . . + n l satisfying relation (1) is called a generalized Γ (L, G) code agreed with the WHM with set of weights v = {r 1 , r 2 , . . . , r l }. l j=1 nj i=1 a (j) i u (j) ′ i (x) u (j) i (x) ≡ 0 mod G(x) , (1)
where deg u

(j) i (x) = r j ≤ deg G(x) , u (j)
i (x) are irreducible polynomials of the degree r j with coefficients from GF (2 m ), (u (j) i (x), G(x)) = 1 for all j = 1, . . . , l; i = 1, . . . , n j , and G(x) is a separable polynomial with coefficients from GF (2 m ).

In order to estimate the maximum power of the set L of the binary generalized Γ (L, G) codes agreed with the WHM, let us give several statements necessary for providing the main result of this paper. Theorem 1. [START_REF] Macwilliams | The Theory of Error-Correcting Codes[END_REF] The number of normalized polynomials of the degree l irreducible over the field GF (2 m ) is determined by the value

I 2 (l) = 1 l d|l µ(d)2 m l d ,
where µ(d) is the Möbius function:

µ(d) =    1 if d = 1 ; (-1) r if d is a product of r different prime numbers ; 0 in all other cases.
It is well-known that the number of unitary separable polynomials of the first degree with coefficients from GF (2 m ) is M 1 2 m = 2 m . Lemma 2. [START_REF] Carlitz | The arithmetic of polynomials in a Galois field[END_REF] The number of unitary separable polynomials of the degree l > 1 with coefficients from the field GF (2 m ) is M l 2 m = 2 ml -2 m(l-1) .

Corollary 1. The number of unitary separable polynomials with coefficients from the field GF (2 m ) whose degrees do not exceed l (l > 1) is equal to

N l 2 m = l i=2 (2 mi -2 m(i-1) ) + 2 m = 2 ml .
The following statement can be obtained from Lemma 1 and Corollary 1 .

Proposition 1. The number of different rational functions

f ′ (x) f (x) , where 1 ≤ deg f (x) ≤ l, f ′ (x) is a formal derivative of f (x) and f (x) is an unitary separable polynomial with coefficients from GF (2 m ), is V l 2 m = 2 ml .

Perfect linear binary codes in the weighted Hamming metric

Let us use definition (4) in order to present a special class of binary generalized Γ (L, G) codes. Now, let l = τ , where τ is a degree of the Goppa polynomial G(x) of the binary generalized Γ (L, G) code. In this case all polynomials of the degree up to τ that are irreducible over GF (2 m ) are chosen as denominators for different elements of the set L. The length of the generalized binary Γ (L, G) code will be equal to n :

n = n 1 + n 2 + . . . + n τ ,
where

n 1 = 2 m , n 2 = I 2 m (2) , . . . , n τ = I 2 m (τ ) -1 .
It means that the first n 1 codeword positions are numerated by all irreducible polynomials of the first degree u (1) i (x). The next n 2 positions are numerated by all irreducible polynomials of the second degree u

(2) i (x) and so on. The last n τ positions are numerated by all irreducible polynomials u (τ ) i (x) of the degree τ with the exception of one that was chosen as Goppa polynomial G(x).

Obviously, when presenting generalized Γ (L, G) codes in the weighted Hamming metric the weight of the i-th code position w i is determined by the degree r i of the polynomial u i (x) used for numbering of this position:

w i = r i . (2)
Let us determine the number of nonzero binary vectors contained in a sphere of radius τ provided that we consider the binary generalized Γ (L, G) code whose position numerators are rational functions

u ′ (j) i (x) u (j) i (x)
, where u (j) i (x)are all the different irreducible polynomials of the degree r (j) ≤ τ with coefficients from GF (2 m ) and G(x) is a irreducible polynomial of degree τ that has not been used for code position numbering. Let e = (e nτ ) be an error vector and E(x) be a syndrome for this vector.

τ j=1 nj i=1 e (j) i u ′ (j) i (x) u (j) i (x) ≡ E(x) mod G(x) , (3)
Lemma 3. Syndromes for various error vectors

e = (e (1) 1 e
(1) 2 . . . e (1) n1 e

(2)

1 . . . e (
u ′ (j) i (x) u (j) i (x) ≡ τ j=1 nj i=1 e (j) i u ′ (j) i (x) u (j) i (x) ≡ E(x) mod G(x) .
is true. This equality can be rewritten in the following form: Hence we can rewrite equation ( 4) in the following form:

τ j=1 nj i=1 (e (j) i ⊕ e (j) i ) u ′ (j) i (x) u (j) i (x) ≡ 0 mod G(x) or ϕ ′ e+b e (x) ϕ e+b e (x) ≡ 0 mod G(x) , (4)
(ψ e+b e (x)) 2 ϕ e+b e (x) ≡ 0 mod G(x) .

To make this equality true, it is necessary to obtain ψ e+b e (x) ≡ 0 mod G(x) , but it is impossible because ψ e+b e (x) = 0 and deg ψ e+b e (x) < τ = deg G(x) .

Since the total number of different nonzero polynomials E(x), deg E(x) < τ with coefficients from GF (2 m ) is equal to 2 mτ -1, and it is equal to the number of different binary vectors e , 0 < wt wH (e) ≤ τ defined by the lemma condition then it follows that for any such vector e there exists an unique polynomial E(x) satisfying equation (3).

Corollary 2. In the weighted Hamming metric, in the sphere of radius τ the number of vectors

e = (e (1) 1 e
(1) 2 . . . e (1) n1 e

. . e (τ ) nτ ), with the length

n = n 1 + n 2 + . . . + n τ ,
where

n 1 = 2 m , n 2 = I 2 m (2) , . . . , n τ = I 2 m (τ ) -1 and v 1 = 1, v 2 = 2, . . . , v τ = τ is equal to W τ n = 2 mτ .
The following statement can be stated by using Corollary 2 and Lemma 1,.

Theorem 2. Binary generalized Γ (L, G) with the Goppa polynomial G(x) (G(x) is an irreducible polynomial with coefficients from GF (2 m ) whose degree is equal to τ ) and the set L where u

(j) i (x)is an irreducible polynomial L = u (j) ′ i (x) u (j)
i (x) i=1,...,nj ;j=1,...,τ of degree j with coefficients from GF (2 m ) and G(x) = u (τ ) i (x) for all i, are the perfect codes in the WHM, i.e., their parameters are on the Hamming bound:

2 k • W τ n = 2 n ,
where n = τ j=1 n j is the codeword length, k is the dimension of the code and W τ n is the number of vectors with length n in a sphere of radius τ in the WHM.

We can obtain the following Corollary from this theorem and Lemma 3. jt j ≤ τ and t j is the number of errors on the n j positions of the codeword.

Remark 1. Note that the binary generalized Γ (L, G) codes perfect in the WHM can be decoded as generalized Γ (L, G) codes. They are decoded by any of the algorithms used for alternant codes, for instance, such as the Euclidian algorithm [START_REF] Macwilliams | The Theory of Error-Correcting Codes[END_REF].

Code example

Let us consider a logical net similar to that given in [START_REF] Deza | Codes Correcting an Arbitrary Sets of Errors[END_REF] whose design is shown in Figure 1. There are n inputs x 1 , x 2 , . . . , x n and 9 outputs y 1 , y 2 , . . . , y 9 . Let us take that only one unit of the net can function erroneously at every time moment [START_REF] Deza | Codes Correcting an Arbitrary Sets of Errors[END_REF]. The faulty operation of unit B 1 leads to an error in output signals y 1 and y 2 so that the corresponding error vector is equal to e 1 = (110000000). An error unit B 2 leads to errors in output signals y 1 and y 3 and the error vector e 2 = (101000000). Errors in other units will lead to similar results. For instance, an error in unit B 7 gives an error in signal y 1 and the error vector will be e 7 = (100000000). Thus, the set of possible error vectors is

F = {e 1 , e 2 , .
. . , e 15 } = {(110000000), (101000000), ( 100100000), (011000000), (010100000), (001100000), ( 100000000), (010000000), (001000000), (000100000), (000010000), (000001000), (000000100), ( 000000010), (000000001)} .

(

) 5
Let us construct a code correcting the total set of error vectors (5) for the logical net. We choose the generalized binary (L, G) code with n 1 = 4, n 2 = 5, k = 5, d wH = 5 that is perfect in the WHM. (It is clear that this code is on the upper bound for linear codes [START_REF] Grassl | Tables of linear Codes and Quantum Codes[END_REF]). Any polynomial of the second degree irreducible over GF (2 2 ) with coefficients from the same field can be chosen as a Goppa polynomial for this code. The total number of such polynomials of the second degree is equal to (2 2 ) 2 -2 2 2 = 6 . One of these polynomials is used as the Goppa polynomial G(x), whereas the rest polynomials form denominators of elements of the set L. Besides, four rational functions with denominators of the first degree The number of codewords of this code is equal to 2 5 . Since the spheres of radius 2 whose centers are in any two codewords of such code do not intersect, then which corresponds to the total number of various binary vectors of length 9.

Table 5 shows all known binary generalized Γ (L, G) codes perfect in the WHM with the length less than 1500.

Notes and comments

In conclusion, it should be noted that the above considered codes are similar to the binary Hamming codes for the HM. The open problem is the existence of perfect codes in the WHM that are similar to non-binary Hamming codes and the Golay codes.

d wH -1 2 ⌋ . Definition 3 .

 23 A binary block code with the length n = n 1 + n 2 + . . . + n l and a set of weights v = (v 1 , v 2 , . . . , v l ) is called perfect in the weighted Hamming metric if its parameters are on the Hamming bound:

τ ) 1 .

 1 . . e (τ ) nτ ), wt wH (e) ≤ τ are different.Proof. We will prove the lemma by contradiction. Assume that we found two different binary vectors e = (e nτ ), wt wH (e) ≤ τ e = ( e

  deg ϕ ′ e+b e (x) < deg ϕ e+b e (x) = wt wH (e ⊕ e) ≤ 2τ. Polynomial ϕ ′ e+b e (x) is the formal derivative in GF (2 m ) and therefore it has nonzero coefficients only at even powers of x. It means that it can be represented as a square of some polynomial ψ e+b e (x): ϕ ′ e+b e (x) = (ψ e+b e (x)) 2 ; deg ψ e+b e (x) < τ .

Corollary 3 .

 3 Binary generalized Γ (L, G) code perfect in the WHM has the minimum distance d ωH = 2τ + 1. It can correct any t ≤ τ j=1 t j errors, where τ j=1

1 x- 1 , 1 x-α , 1 x-α 2 , 1 x 2 Table 1 .

 1111121 are used in the set L, where α is a primary element of GF (2 2 ). Let us find the number of vectors in the sphere of radius 2n1 = 4 n2 = 5 n = 9 v1 = 1 v2 = 2 t1 t2 t = t1 + t2 0 Binary generalized Γ (L, G) code with (n1 = 4, n2 = 5, k = 5, dwH = 5)in the weighted Hamming metric for this code:

Fig. 1 .

 1 Fig. 1. Logical net


Adrien Bernhardt

André Maximo

Luiz Velho

Houssam Hnaidi

Marie-Paule Cani

Real-time Terrain Modeling using CPU-GPU Coupled Computation

Introduction and Related Work

In many editing tools, especially sketch-based modeling, it is important to have real-time feedto help improve the editing quality. This importance is emphasized particularly in sketch-based terrain modeling, being able to see the terrain morphing at the same time the drawing edition occurs constitutes a great user experience (see Figure 1(a)). In this work, we propose a real-time terrain modeling tool by combining a fast GPU-based terrain solver [START_REF] Hnaidi | Feature based terrain generation using diffusion equation[END_REF]] with a lightweight CPU-based data structure. Our tool is capable of dynamically generate multi-resolution heightmaps, enabling it to tessellate different parts of the terrain at different resolutions.

In our framework, we have two types of editing interactions: in the spirit of [START_REF] Gain | Terrain Sketching[END_REF]], the user can draw strokes creating elevations and crevices; and previous strokes can be interactively moved to different regions of the terrain. Differently than Gain et al.'s system, we do not extract noise from the user strokes to make the terrain more realistic, however we use a CPU-GPU coupled method to drastically improve the performance of our tool, generating terrains two orders of magnitude faster than Gain et al.'s work.

The terrain modeling in our approach is accomplished by combining the multi-grid GPU terrain solver of [START_REF] Hnaidi | Feature based terrain generation using diffusion equation[END_REF]] with an adaptive tessellation-based rendering shader capable of handling dynamic heightmaps. The main contribution of Hnaidi et al.'s work is to propose a GPU-based multi-grid diffusion equation solver, which interpolates not only heights but also amplitude and frequency noise. Our modeling tool uses Hnaidi et al.'s solver to allow an interactive manipulation of complex terrain primitives.

Real-time Terrain Modeling

Creation of terrain models in realtime involves dealing with dynamically changing data that increases exponentially depending on the terrain resolution. In order to provide a real-time terrain modeling tool, we make use of two complementary approaches. First, a coarse version of the terrain is maintained in the CPU using a quadtree (see Figure 1(b)), where regions closer to the viewer are subdivided more than far regions. This simple and lightweight data structure fits the CPU main role of data control, while allowing it to send adaptive quad primitives to the GPU. Second, a fine version of the terrain is produced in the GPU using the tessellation control and evaluation shaders. The tessellation control shader is responsible to subdivide regularly each patch primitive, i.e. the quad leaf node sent by the CPU, while the tessellation evaluation shader reads the height values from a texture. The first step of our algorithm is to update the quadtree data structure using a LOD-based approach. We consider the projection of the bounding box of each quad node by reading the minimum and maximum height value that falls inside the node. The projection is used to determine if the patch node needs to be sent to the GPU and to interactively adapt the quadtree in a way that each projection shape has about the same size.

The second step of our algorithm is to translate the sketch-defined terrain primitives to constraints that are used by the multi-grid GPU solver of [START_REF] Hnaidi | Feature based terrain generation using diffusion equation[END_REF]. This solver provides a sequence of increasing resolution grids, up to an arbitrary size, which we store in a mipmap-pyramid texture to be read in our tessellation control and evaluation shaders. The tessellation control shader uses it to decide the subdivision level of each patch, while the tessellation evaluation shader uses it to place each generated vertex at the proper height value. The multi-resolution texture is also used by the CPU, but only a small resolution (6th mipmap level) of it since the quadtree minimum leaf size is much bigger than the texel from the highest resolution texture.

In our early experiments we discovered that by modifying several aspects of the original GPU solver, we are able to compute the entire multi-resolution heightmap texture (with 4K × 4K maximum size) in 80 ms using an off-the-shelf graphics card. Moreover, we use the CPU to control the solver iterations and stop at a certain resolution and then resume computing when the GPU is idle. Another interesting feature of our method is the balance between terrain generation in the CPU and in the GPU, we can control this balance by simply changing the quadtree refinement. With these features, the user can draw strokes and see at the same time the terrain morphing to the drawing. Terrain primitives, such as cliffs and mountains, are controlled seamlessly in our framework. In short, we believe this work has the potential to become an effective terrain modeling tool, creating high-quality terrain models in real-time.

Figure 1 :

 1 Figure 1: Terrain modeling example: as each stroke is drawn (a) or manipulated (c), the terrain is tessellated in the GPU to follow the stroke. In the CPU, the quadtree data structure (b) controls the quad patches sent to the GPU.