Saverio Morfu
email: smorfu@u-bourgogne.fr

Patrick Marquié

Brice Nofiele

Dominique Ginhac

B Nofiélé

Nonlinear Systems for Image Processing

published or not. The documents may come   L'archive ouverte pluridisciplinaire

I. INTRODUCTION

Since about one hundred years, nonlinear science has attracted the attention of researchers to circumvent the limitation of linear theories in the explanation of natural phenomenons. Indeed, nonlinear differential equations can model the behavior of oceans surface [START_REF] Scott | Nonlinear Science emergence and dynamics of coherent structures[END_REF], the recurrence of ice age [START_REF] Benzi | Stochasic resonance in climatic change[END_REF], the transport mechanisms in living cells [START_REF] Murray | Mathematical Biology[END_REF], the information transmission in neural networks [START_REF] Nagumo | An active pulse transmission line simulating nerve axon[END_REF][START_REF] Scott | Nonlinear Science emergence and dynamics of coherent structures[END_REF]Izhikevich, 2007), the blood pressure propagation in arteries (Paquerot and Remoissenet, 1994) or the excitability of cardiac tissues [START_REF] Beeler | Reconstruction of the action potentials of ventricular myocardial fibers[END_REF][START_REF] Keener | Propagation and its failure in coupled systems of discrete excitable cells[END_REF]. Therefore, nonlinear science appears as the most important frontier for a better understanding of nature [START_REF] Remoissenet | Waves called solitons: concepts and experiments[END_REF]. In the recent field of engineering science [START_REF] Zakharov | Optical Solitons: Theoretical Challenges and Industrial perspectives[END_REF][START_REF] Agrawal | Fiber-Optic Communication systems[END_REF], taking into account the nonlinearity has allowed to achieve spectacular progresses in terms of transmission capacities in optical fibers via the concept of soliton [START_REF] Remoissenet | Waves called solitons: concepts and experiments[END_REF]. More recently, nonlinear differential equation, arising in many areas of physics, biology, chemistry and ecology, have naturally inspired unconventional methods of processing which allow to transcend the limitations of classical linear methods [START_REF] Teuscher | Proceedings of the 2005 Workshop on Unconventional Computing, From cellular Automata to Wetware[END_REF]. This growing interest for processing applications based on the properties of nonlinear systems can be explained by the observation that fundamental progress in several fields of computer science seems sometimes to stagnate. Novel ideas coming from interdisciplinary fields often open new directions of research with unsuspected applications [START_REF] Teuscher | Proceedings of the 2005 Workshop on Unconventional Computing, From cellular Automata to Wetware[END_REF].

On the other hand, complex processing tasks require intelligent systems that are able to adapt and learn by mimicking the behaviour of human brain. Biologically inspired systems, most often described by nonlinear reaction-diffusion equations, have then been proposed as convenient solutions to solve very complicated problems unaccessible to modern von Neumann computers. It was in this context that the concept of CN N has been introduced by L. Chua and L. Yang as a novel class of information processing systems with potential applications in such areas as image processing and pattern recognition (Chua and Yang, 1988;Chua and Yang, 1988). In fact, CN N is an acronym for Cellular Neural Network when used in the context of Brain science or Cellular Nonlinear Network in the context of emergence and complexity [START_REF] Chua | CNN: a paradigm for complexity[END_REF]. Since the pioneer work of L.O Chua, the CNN paradigm has rapidly evolved to cover a wide range of applications drawn from numerous disciplines, including artificial life, biology, chemistry, physics, information science, nonconventional methods of computing [START_REF] Holden | Can excitable media be considered as computational systems?[END_REF], video coding [START_REF] Venetianer | Analogic CNN Algorithms for some Image compression and restoration tasks[END_REF][START_REF] Arena | An object oriented segmentation on analog CNN chip[END_REF], quality control by visual inspection [START_REF] Occhipinti | Defects detection and characterization by using cellular neural networks[END_REF], cryptography [START_REF] Yu | Cryptography based on delayed chaotic neural networks[END_REF][START_REF] Caponetto | Sc-CNNs for chaotic signal applications in secure communication systems[END_REF], signal-image processing [START_REF] Julián | A piecewise-linear simplicial coupling cell for CNN Gray-Level Image Processing[END_REF],... (see [START_REF] Tetzlaff | Cellular Neural Networks and their applications[END_REF] for an overview of the applications).

In summary, these latest two decades devoted to the study of CN N s have led scientists to solve problems of artificial intelligence by combining the highly parallel multiprocessor architecture of CN N s with the properties inherited from the nonlinear bio-inspired systems. Among the tasks of high computational complexity routinely performs with nonlinear systems, one can cite finding the optimal path in a two dimensional vector field [START_REF] Agladze | Finding the optimal path with the aid of chemical wave[END_REF], image skeletonisation [START_REF] Chua | CNN: a paradigm for complexity[END_REF], finding the shortest path in a labyrinth [START_REF] Chua | CNN: a paradigm for complexity[END_REF][START_REF] Rambidi | Chemical reaction-diffusion implementation of finding the shortest paths in a labyrinth[END_REF], or controlling mobile robot [START_REF] Adamatzky | Experimental implementation of mobile robot taxis with onboard Belousov-Zhabotinsky chemical medium[END_REF]. However, the efficiency of these nonlinear systems for signal-image processing or pattern recognition does not come only from their biological background. Indeed, the nonlinearity offers an additional dimension lying in the signal amplitude, which gives rise to novel properties not shared by linear systems. Noise removal with a nonlinear dissipative lattice [START_REF] Marquié | Diffusion effects in a nonlinear electrical lattice[END_REF][START_REF] Comte | Noise removal using a twodimensional diffusion network[END_REF], contrast enhancement based on nonlinear oscillators properties [START_REF] Morfu | A nonlinear oscillators network devoted to image processing[END_REF], edge detection exploiting vibration noise (Hongler et al, 1998), optimization by noise of non-optimum problems or signal detection aided by noise via the famous stochastic resonance phenomenon [START_REF] Gammaitoni | Stochastic Resonance[END_REF][START_REF] Chapeau-Blondeau | Stochastic resonance and the benefit of noise in nonlinear systems, Noise, Oscillators and Algebraic Randomness-From Noise in Communication Systems to Number Theory[END_REF][START_REF] Comte | Stochastic Resonance: another way to retrieve subthreshold digital data[END_REF], constitute a non restrictive list of spectacular examples where the properties of nonlinear systems have allowed to overcome the limitation of classical linear approaches.

Owing to the rich variety of potential applications inspired by nonlinear systems, rapidly the efforts of researchers have focused on the experimental realization of such efficient information processing devices. Two different strategies were introduced [START_REF] Kuhnert | A new optical photochemical device in a light-sensitive chemical active medium[END_REF]Chua and Yang, 1988) and nowadays, the fascinating challenge of artificial intelligence implementation with CN N is still under investigation.

The first technique dates back from the late eighties with the works of L. Kuhnert who proposed to take benefit of the properties of Belousov-Zhabotinsky type media for image processing purposes [START_REF] Kuhnert | A new optical photochemical device in a light-sensitive chemical active medium[END_REF][START_REF] Kuhnert | Image processing using light-sensitive chemical waves[END_REF]. The main idea is that each mico-volume of the active photo-sensitive chemical medium acts as a one-bit-processor corresponding to reduced/oxidized state of the catalyst [START_REF] Agladze | Finding the optimal path with the aid of chemical wave[END_REF]. This feature of chemical photosensitive nonlinear media has allowed to implement numerous tools of image processing. Edge enhancement, classical operations of mathematical morphology, the restoration of individual components of an image with overlapped components [START_REF] Rambidi | Image processing using light-sensitive chemical waves[END_REF], the image skeletonisation [START_REF] Adamatzky | Experimental reactiondiffusion pre-processor for shape recognition[END_REF], the detection of urban roads or the analysis of medical images [START_REF] Teuscher | Proceedings of the 2005 Workshop on Unconventional Computing, From cellular Automata to Wetware[END_REF] represent a brief overview of processing tasks computed by chemical nonlinear media. However, even considering the large number of chemical "processors", the very low velocity of trigger waves in chemical media is sometimes incompatible with real time processing constraints imposed by practical applications [START_REF] Agladze | Finding the optimal path with the aid of chemical wave[END_REF]. Nevertheless, the limitations of these unconventional methods of computing no way dismiss the efficiency and high prospects of the processing developed with active chemical media [START_REF] Adamatzky | On some limitations of reaction-diffusion chemical computers in relation to Voronoi diagram and its inversion[END_REF].

By contrast, analog circuits do not share the weakness of the previous strategy of integration. Therefore, because of their real-time processing capability, electronic hardware devices constitute the most common way to implement CN N s (Chua and Yang, 1988). The first step to electronically develop a CN N for image processing purposes consists of designing an elementary cell. More precisely, this basic unit of CN N s usually contains linear capacitors, linear resistors, linear and nonlinear controlled sources (Chua and Yang, 1988;Comte and Marquié, 2003). Next, to complete the description of the network, a coupling law between cells is introduced. Owing to the propagation mechanism inherited from the continuous-time dynamics of the network, the cells do not only interact with their nearest neighbors but also with cells that are not directly connected together. Among the applications which can be electronically realized, one can cite character recognition (Chua and Yang, 1988), edge filtering (Comte et al, 2001;[START_REF] Chen | Image processing algorithms realized by discrete-time cellular neural networks and their circuit implementation[END_REF], noise filtering [START_REF] Marquié | Diffusion effects in a nonlinear electrical lattice[END_REF][START_REF] Comte | Noise removal using a twodimensional diffusion network[END_REF][START_REF] Julián | A piecewise-linear simplicial coupling cell for CNN Gray-Level Image Processing[END_REF], contrast enhancement and gray level extraction with a nonlinear oscillators network [START_REF] Morfu | Image processing with a cellular nonlinear network[END_REF][START_REF] Morfu | Experiments on an electrical nonlinear oscillators network[END_REF]. On the other hand, the principle of CN N s integration with discrete electronic components is closely related to the development of nonlinear electrical transmission lines (NLTLs) [START_REF] Remoissenet | Waves called solitons: concepts and experiments[END_REF]. Indeed, under certain conditions [START_REF] Chua | CNN: a paradigm for complexity[END_REF], the parallel processing of information can be ruled by nonlinear differential equations which also describe the evolution of the voltage at the nodes of an electrical lattice. It is then clear that considering a one dimensional lattice allows signal filtering, while extending the concept to a two dimensional network can provide image processing applications.

The development of NLTLs was mainly motivated by the fact that these systems are quite simple and relatively un expansive experimental devices allowing to study quantitatively the properties of nonlinear waves [START_REF] Scott | Active and Nonlinear Wave Propagation in Electronics[END_REF]. In particular, since the pioneering works by Hirota and Suzuki [START_REF] Hirota | Studies on lattice solitons by using electrical networks[END_REF] and by Nagashima and Amagishi [START_REF] Nagashima | Experiments on the Toda lattice using nonlinear transmission line[END_REF] on electrical lines simulating the Toda lattice [START_REF] Toda | Wave propagation in anharmonic lattices[END_REF], these NLTLs, which can be considered as analog simulators, provide a useful way to check how the excitations behave inside the nonlinear medium [START_REF] Jäger | Characteristics of travelling waves along the nonlinear transmission line for monolithic integrated circuits: a review[END_REF][START_REF] Kuusela | Soliton experiments in transmission lines[END_REF][START_REF] Marquié | Observation of nonlinear localized modes in an electrical lattice[END_REF][START_REF] Yamgoué | Noise effects on gap wave propagation in a nonlinear discrete LC transmission line[END_REF].

This chapter is mainly devoted to the presentation of some particular nonlinear processing tools and to discuss their electronic implementation with discrete components.

After having introduced the readers into a brief mechanical description of nonlinear systems, we first present a review of the properties of both purely inertial systems and overdamped systems. Then, in the following sections, we take advantage of these properties to develop unconventional methods of processing. Especially, considering the features of purely inertial systems, we reveal the possibility to perform various image processing tasks such as contrast enhancement of a weakly contrasted picture, the extraction of gray levels, or the encryption of an image. The electronic sketch of the elementary cell of this inertial CN N is proposed and the nonlinear properties which allows the previous image processing tasks are experimentally investigated. Next, the third part of this chapter is exclusively devoted to the filtering applications inspired by reaction-diffusion media, like for instance, noise filtering, edge detection or extraction of interest regions in a weakly noisy contrasted picture. In each case, the elementary cell of the electronic CN N is developed and we experimentally investigate its behavior in the specific context of signal-image processing. We conclude by discussing the possible microelectronic implementations of the previous nonlinear systems. In addition, the last section contains some perspectives for future developments inspired by recent properties of nonlinear systems. In particular, we present a paradoxical nonlinear effect known as stochastic resonance [START_REF] Benzi | Stochasic resonance in climatic change[END_REF][START_REF] Gammaitoni | Stochastic Resonance[END_REF][START_REF] Chapeau-Blondeau | Noise assisted propagation over a nonlinear line of threshold elements[END_REF] which is supposed to have potential applications in visual perception [START_REF] Simonotto | Visual perception of stochastic resonance[END_REF].

We trust that the multiple topics proposed in this contributions will help the readers in better understanding the potential applications based on the properties of nonlinear systems. Moreover, the various electronic realizations presented in this chapter will constitute a serious background for future experiments and studies devoted to nonlinear phenomena. Written for an interdisciplinary readership of physicist and engineers, we finally hope that this chapter will encourage the readers to perform their own experiments.

II. MECHANICAL ANALOGY

In order to understand the image processing tools inspired by the properties of nonlinear systems, we present here a mechanical analogy of these nonlinear systems. From a mechanical point of view, we consider a chain of particles of mass M submitted to a nonlinear force f deriving from a potential Φ and coupled with springs of strength D. If W n represents the displacement of the particle n, the fundamental principle of the mechanics writes:

M d 2 W n dt 2 + λ dW n dt = - dΦ dW n + R n , (1)
where M d 2 W dt 2 represents the inertia term and λ dW dt corresponds to a friction force. Furthermore, the resulting elastic force R n applied to the n th particle by its neighbors can be defined by:

R n = D j∈N r W j -W n , (2)
where N r is the neighborhood, namely N r = {n -1, n + 1} in the case of a one dimensional chain.

We propose to investigate separately the purely inertial case, that is

M d 2 W dt 2 >> λ dW dt
, and the overdamped one deduced when M d 2 W dt 2 << λ dW dt .

A. Overdamped Case

In this section, an overdamped system is presented by neglecting the inertia term of eq. ( 1) compared to the friction force. We specifically consider λ = 1 and the case of a cubic nonlinear force

f (W ) = -W (W -α)(W -1), (3)
deriving from the double well potential Φ(W ) = -W 0 f (u)du as represented in figure 1 for different worths of α. The roots of the nonlinear force 0, and 1 correspond to the positions of the local minima of the potential, namely the well bottoms, whereas the root α represents the position of the potential maximum. The nonlinearity threshold α defines the potential barrier ∆ between the potential minimum with the highest energy and the

F ( W ) (b)
FIGURE 1. Double well potential deduced from the nonlinear force (3). (a) For α < 1/2 the well bottom with highest energy is located at W = 0, the potential barrier is given by

∆ = α 0 f (u)du = φ(α) -φ(0) . (b) For α > 1/
2 the symmetry of the potential is reversed: W = 1 becomes the position of the well bottom of highest energy and the potential barrier is

∆ = α 1 f (u)du = φ(α) -φ(1).
potential maximum. To explain the propagation mechanism in this chain, it is convenient to define the excited state by the position of the potential minimum with the highest energy, and the rest state by the position corresponding to the minimum of the potential energy. As shown in figure 1.(a), the excited state is 0 and the rest state is 1 when the nonlinearity threshold α < 1/2. In the case α > 1/2, since the potential symmetry is reversed, the excited state becomes 1 and the rest state is 0 (figure 1.(b)). The equation which rules this overdamped nonlinear systems can be deduced from eq. ( 1). Indeed, when the second derivative versus time is neglected compared to the first derivative and when λ = 1, eq. ( 1) reduces to the discrete version of Fisher's equation, introduced in the 1930's as a model for genetic diffusion [START_REF] Fisher | The wave of advance of advantageous genes[END_REF]:

dW n dt = D(W n+1 + W n-1 -2W n ) + f (W n ). (4)

Uncoupled case

We first investigate the uncoupled case, that is D = 0 in eq. ( 4) to reveal the bistability of the system. The behavior of a single particle of displacement W and initial position W 0 obeys to

dW dt = -W (W -α)(W -1). ( 5
)
The zeros of the nonlinear force f , W = 1 and W = 0 correspond to stable steady states, whereas the state W = α is unstable. The stability analysis can be realized by solving eq. ( 5) substituting the nonlinear force f = -W (W -α)(W -1) by its linearized expression near the considered steady states W * ∈ {0, 1, α}. If f W (W * ) denotes the derivative versus W of the nonlinear force for W = W * , we are led to solve:

dW dt = f W (W * )(W -W * ) + f (W * ). (6)
The solution of Eq. ( 6) can then be easily expressed as

W (t) = W * + Ce fW (W * )t - f (W * ) f W (W * ) , ( 7
)
where C is a constant depending on the initial condition, that is the initial position of the particle. The solution [START_REF] Yu | Cryptography based on delayed chaotic neural networks[END_REF], obtained with a linear approximation of the nonlinear force f , shows that the stability is set by the sign of the argument of the exponential function. Indeed, for W * = 0 and W * = 1, the sign of f W (W * ) is negative, involving that W (t → ∞) tends to a constant. Therefore, the two points W * = 0 and W * = 1 are stable steady states.

On the other hand, for

W * = α, f W (W * ) is positive, inducing a divergence for W (t → ∞). W * = α is an unstable steady state.
We now focus our attention to the particular case α = 1/2 since it will allow interesting applications in signal and image processing context.

This case is intensively developed in Appendix A, where it is shown that the displacement of a particle with initial position W 0 can be expressed by

W (t) = 1 2 1 + W 0 -1 2 (W 0 -1 2 ) 2 -W 0 (W 0 -1)e -1 2 t . ( 8
)
This theoretical expression is compared in figure 2 to the numerical results obtained solving eq. ( 5) using a fourth order Runge-Kutta algorithm with integrating time step dt = 10 -3 . As shown in figure 2, when the initial condition W 0 is below the unstable state α = 1/2, the particle evolves toward the steady states 0. Else, if the initial condition W 0 exceeds the unstable state α = 1/2, the particle evolves towards the other steady state 1. Therefore, the unstable states α = 1/2 acts as a threshold and the system exhibits a bistable behaviour.

Coupled case

We now consider the coupled case (D = 0). In such systems ruled by eq. ( 4), the balance between the dissipation and the nonlinearity gives rise to the propagation of a kink, that is a localized wave, called diffusive soliton which propagates with constant velocity and profile [START_REF] Remoissenet | Waves called solitons: concepts and experiments[END_REF]. To understand the propagation mechanism, we first consider the weak coupling limit and the case α < 1/2. The case of strong couplings, which corresponds to a continuous medium, will be discussed later since it allows to theoretically characterize the waves propagating in the medium.

1.1. Weak coupling limit As shown in figure 3.(a), initially, all particles of the chain are located at the position 0, that is at the excited state. To initiate a kink, an external forcing allows the first particle to cross the potential barrier in W = α and to fall in the right well, at the rest state defined by the position W = 1. Thanks to the spring coupling the first particle to the second one, but despite the second spring, the second particle attempts to cross the potential barrier with height

∆(α) = - α 4 12 + α 3 6 (Morfu, 2003) (see figure 3.(b)
). According to the value of the resulting force applied to the second particle by the two springs compared to the nonlinear force f between [0, α[, two behaviors may occur:

1. If the resulting elastic force is sufficiently important to allow the second particle to cross the potential barrier ∆(α), then this particle fall in the right well and pulls the next particle down in its fall. Since each particle

D ( a ) 0 1 W n W 4 W 3 W 2 W 1 D D D f ( W n ) ( t = 0 ) D a (a) D ( a ) 0 W 2 a 1 W n W 4 W 3 W 2 W 1 D D D f ( W n ) ( t > 0 ) ? (b) FIGURE 3. Propagation mechanism. (a)
Initially all particles of the chain are in the excited state 0, that is at the bottom of the well with highest energy. (b) State of the chain for t > 0. The first particle has crossed the potential barrier ∆ and attempts to pull the second particle down in its fall.

of the chain successively undergoes a transition from the excited state 0 to the rest state 1, a kink propagates in the medium. Moreover, its velocity increases versus the coupling and as the barrier decreases (namely, as α decreases). 2. Else, if the resulting force does not exceed a critical value,(i.e. if D < D * (α)), the second particle cannot cross the potential barrier and thus stays pinned at a position w in [0; α[: it is the well known propagation failure effect [START_REF] Keener | Propagation and its failure in coupled systems of discrete excitable cells[END_REF][START_REF] Erneux | Propagating waves in discrete bistable reactiondiffusion systems[END_REF][START_REF] Kladko | Universal scaling of wave propagation failure in arrays of coupled nonlinear cells[END_REF]Comte et al, 2001).

The mechanical model associated with eq. ( 4) reveals that in the weak coupling limit the characteristics of the nonlinear system are ruled by the coupling D and the nonlinear threshold α. Moreover, the propagation of a kink is due to the transition from the excited state to the rest state and is only possible when the coupling D exceeds a critical value D * (α).

Limit of continuous media

The velocity of the kink and its profile can be theoretically obtained in the limit of continuous media, that is when the coupling D is large enough compared to the nonlinear strength.

Then, in the continuous limit, the discrete Laplacian of eq. ( 4) can be replaced by a second derivative versus the space variable z:

∂W ∂t = D ∂ 2 W ∂z 2 + f (W ). (9)
This equation, introduced by Nagumo in the 1940s as an elementary representation of the conduction along an active nerve fiber, has an important meaning in understanding transport mechanism in biological systems [START_REF] Nagumo | An active pulse transmission line simulating nerve axon[END_REF][START_REF] Murray | Mathematical Biology[END_REF]. Unlike the discrete equation ( 4), the continuous equation ( 9) admits propagative kink solution only if 1 0 f (u)du = 0, which reduces to α = 1/2 in the case of the cubic force (3) [START_REF] Scott | Nonlinear Science emergence and dynamics of coherent structures[END_REF]. Introducing the propagative variable ξ = z -ct, these kinks and anti-kinks have the form [START_REF] Fife | Mathematical Aspects of Reacting and Diffusing Systems[END_REF][START_REF] Henry | Geometric theory of semilinear parabolic equations[END_REF])

W (ξ) = 1 2 1 ± tanh 1 2 √ 2D (ξ -ξ 0 ) , (10)
where ξ 0 is the initial position of the kink for t = 0 and where the kink velocity is defined by c = ± D/2(1 -2α).

When α < 1/2, the excited state is 0, and the rest state is 1. Therefore, the rest state 1 spreads in the chain, which set the sign of the velocity according to the profile of the kink initiated in the nonlinear system:

1. If the profile is given by W When α > 1/2, since the symmetry of the potential is reversed , the excited states becomes 1 and the rest state is 0. The propagation is then due to a transition between 1 to 0 which provides the following behaviour:

(ξ) = 1 2 1 -tanh 1 2 √ 2D (ξ -ξ 0 ) ,
1. If W (ξ) = 1 2 1 -tanh 1 2 √ 2D (ξ -ξ 0 ) ,

Else if

W (ξ) = 1 2 1 + tanh 1 2 √ 2D (ξ -ξ 0 ) ,

B. Inertial systems

In this section, we neglect the dissipative term of eq. ( 1) compared to the inertia term and we restrict our study to the uncoupled case. Moreover, in image processing context, it is convenient to introduce a nonlinear force f under the form

f (W ) = -ω 2 0 (W -m)(W -m -α)(W -m + α), (11)
where, m and α < m are two parameters which allow to adjust the width and height ∆ = ω 2 0 α 4 /4 of the potential Φ (fig. 5): The nonlinear differential equation which rules the uncoupled chain can be deduced by inserting the nonlinear force (11) into eq. ( 1) with D = 0. Neglecting the dissipative term, the particles of unitary mass are then ruled by the following nonlinear oscillator equations:

Φ(W ) = - W 0 f (u)du. (12)
d 2 W i dt 2 = f (W i ). (13)

Theoretical analysis

We propose here to determine analytically the dynamics of the nonlinear oscillators obeying to eq. ( 13) [START_REF] Morfu | A nonlinear oscillators network devoted to image processing[END_REF][START_REF] Morfu | Contrast enhancement with a nonlinear oscillators network[END_REF]. Setting x i = W i -m, eq. ( 13) can be rewritten as

d 2 x i dt 2 = -ω 2 0 x i (x i -α)(x i + α). (14)
Noting x 0 i the initial position of the particle i and considering that all the particles have initially a null velocity, the solutions of eq. ( 14) can be expressed with the Jacobian elliptic functions as

x i (t) = x 0 i cn(ω i t, k i ), (15)
where ω i and 0 ≤ k i ≤ 1 represent respectively the pulsation and the modulus of the cn function (see recall on the properties of Jacobian elliptic function in Appendix B). Deriving twice eq. ( 15) and using the properties (85), we get

dx i dt = -x 0 i ω i sn(ω i t, k i )dn(ω i t, k i ), d 2 x i dt 2 = -x 0 i ω 2 i cn(ω i t, k i ) dn 2 (ω i t, k i ) -k i sn 2 (ω i t, k i ) . ( 16
)
Using the identities ( 86) and (87), eq. ( 16) can be rewritten as

d 2 x i dt 2 = - 2k i ω 2 i x 0 2 i x x 2 - 2k i -1 2k i x 0 2 i . ( 17
)
Identifying this last expression with eq. ( 14), we get the pulsation of the Jacobian elliptic function

ω i = ω 0 x 0 2 i -α 2 , (18)
and its modulus

k i = 1 2 x 0 2 i x 0 2 i -α 2 . ( 19
)
Lastly, introducing the initial condition W 0 i = x 0 i + m, the solution of eq. ( 13) can be straightforwardly deduced from eqs. (15), ( 18) and ( 19):

W i (t) = m + (W 0 i -m)cn(ω i t, k i ), (20)
with

ω i (W 0 i ) = ω 0 (W 0 i -m) 2 -α 2 and k i (W 0 i ) = 1 2 (W 0 i -m) 2 (W 0 i -m) 2 -α 2 . ( 21
)
Both the modulus and the pulsation are driven by the initial condition W 0 i . Moreover, the constraints to ensure the existence of the pulsation ω i and of the modulus respectively writes (W 0 i -m) 2 -α 2 ≥ 0 and 0 ≤ k i ≤ 1. These two conditions restrict the range of the allowed initial conditions the pulsation and the modulus are represented versus the initial condition W 0 i . Note that this allowed range of initial conditions corresponds also to a particle with an initial potential energy exceeding the barrier ∆ between the potential extrema (see figure 5).

W 0 i to -∞; m -α √ 2 m + α √ 2; + ∞ ,

Nonlinear oscillators properties

To illustrate the properties of nonlinear oscillators, we consider a chain of length N = 2 particles with a weak difference of initial conditions and with a null initial velocity. The dynamics of these two oscillators is ruled by eq. ( 20), where the pulsation and modulus of both oscillators are driven by their respective initial condition. Moreover, we have restricted our study to the case of the following nonlinearity parameters m = 2.58, α = 1.02, ω 0 = 10 4 . We have applied the initial condition W 0 1 = 0 to the first oscillator, while the initial condition of the second oscillator is set to W 0 2 = 0.2, which corresponds to the situation of fig. 5. Figure 7.(a) shows that the oscillations of both particles take place in the

0 0 . 5 1 1 . 5 2 2 . 5 3 0 1 2 3 4 5 6 t i m e W 1 ( t ) 0 1 2 3 4 5 W 2 ( t ) 0 0 . 5 1 1 . 5 2 2 . 5 3 t i m e t o p t
x ( 1 0 -3 )

x ( 1 0 -3 ) (a) range [W 0 i ; 2m -W 0 i ] as predicted by eq. ( 20), that is [0; 5.16] for the first oscillator and [0; 4.96] for the second one. Moreover, owing to their difference of initial amplitude and to the nonlinear behavior of the system, the two oscillators quickly attain a phase opposition for the first time at t = t opt = 1.64 × 10 -3 . This phase opposition corresponds to the situation where the first oscillator has reached its minimum W 1 (t opt ) = 0, whereas the second oscillator has attained its maximum W 2 (t opt ) = 4.96. As shown in figure 7.(b), the displacement difference δ(t) = W 2 (t) -W 1 (t) is then maximum for t = t opt and becomes δ(t opt ) = 4.96. For this optimal time, a "contrast enhancement" of the weak difference of initial conditions is realized, since initially the displacement difference was δ(t = 0) = 0.2. Note that in fig. 7.(b), the displacement difference between the two oscillators also presents a periodic behavior with local minima and local maxima. In particular, the difference δ(t) is null for t = 3.96 × 10 -5 , t = 1.81 × 10 -4 , t = 3.5 × 10 -4 , t = 5.21 × 10 -4 ; minimum for t = 1.4 × 10 -4 , t = 4.64 × 10 -4 , t = 1.47 × 10 -3 and maximum for t = 3 × 10 -4 , t = 6.29 × 10 -4 , t = 1.64 × 10 -3 . These characteristic times will be of crucial interest in image processing context to define the filtered tasks performed by the nonlinear oscillators network. Figure 6.(a) reveals that the maximum variation of the pulsation compared to the amplitude W 0 i , that is ∆ω/ω 0 , is reached for

-5 -4 -3 -2 -1 0 1 2 3 4 5 0 0 . 5 1 1 . 5 2 2 . 5 3 t i m e d ( t ) = W 2 ( t ) -W 1 ( t ) t o p t x ( 1 0 -3 ) (b)
W 0 i = m -α √ 2
, that is for a particle with an initial potential energy near the barrier ∆. Therefore, to quickly realize a great amplitude contrast between the two oscillators, it could be interesting to launch them with an initial amplitude near m -α √ 2, or to increase the potential barrier height ∆. We choose to investigate this latter solution by tuning the parameter of the nonlinearity α, when the initial amplitude of both oscillators remains W 0 1 = 0 and W 0 2 = 0.2. The results are reported in figure 8, where we present the evolution of the difference δ(t) for different worths of α.

0 1 2 3 x 1 0 -3 -5 -3 -1 1 3 5 t d ( t ) = W 2 ( t ) -W 1 ( t ) t o p t (a) t o p t 0 1 2 3 x 1 0 -3 -5 -3 -1 1 3 5 t d ( t ) = W 2 ( t ) -W 1 ( t ) (b) 0 1 2 3 x 1 0 -3 -5 -3 -1 1 3 5 t d ( t ) = W 2 ( t ) -W 1 ( t ) t o p t (c) t o p t -5 -3 -1 1 3 5 0 1 2 3 x 1 0 -3 t d ( t ) = W 2 ( t ) -W 1 ( t ) (d)
As expected, when the nonlinearity parameter α increases, the optimal time is significatively reduced. However, when α is adjusted near the critical value (m -W 0 2 )/ √ 2 as in figure 8.(d), the optimum reached by the difference δ(t) is reduced to 4.517 for α = 1.63 instead of 4.96 for α = 1.02. Even if it is not the best contrast enhancement which can be performed by the system, the weak difference of initial conditions between the two oscillators is nevertheless strongly enhanced for α = 1.63.

To highlight the efficiency of nonlinear systems, let us consider the case of a linear force f (W ) = -ω 0 W in eq. ( 13). In the linear case, the displacement difference δ(t) between two harmonic oscillators can be straightforwardly expressed as

δ(t) = ǫ cos(ω 0 t), (22)
where ǫ represents the slight difference of initial conditions between the oscillators. This last expression shows that it is impossible to increase the weak difference of initial conditions since the difference δ(t) always remains in the range [-ǫ; ǫ]. Therefore, taking into account nonlinearity is a convenient solution to overcome the limitation of linear system and to enhance a weak amplitude contrast.

III. INERTIAL SYSTEMS

In this section, we present different image processing tasks inspired by the properties of the nonlinear oscillators presented in section II.B. Their electronic implementation is also discussed.

A. Image processing

By analogy with a particle experiencing a double well potential, the pixel number (i, j) is analog to a particle (oscillator) whose initial position corresponds to the initial gray level W 0 i,j of this pixel. Therefore, if N ×M denotes the image size, we are led to consider a two dimensional network, or Cellular Nonlinear Network (CN N ), consisting of uncoupled nonlinear oscillators. The node i, j of this CNN obeys to

d 2 W i,j dt 2 = -ω 2 0 (W i,j -m -α)(W i,j -m + α)(W i,j -m), (23)
with i = 1, 2...N and j = 1, 2.., M . Note that we take into account the range of oscillations [0; 2m -W 0 i,j ] predicted in section II.B.2 to define the gray scale of the images, namely 0 for the black level and 2m = 5.16 for the white level.

The image to process is first loaded as initial condition at the nodes of the CN N . Next, the filtered image for a processing time t can be deduced noting the position reached by all oscillators of the network at this specific time t. More precisely, the state of the network at a processing time t is obtained by solving numerically eq. ( 23) with a fourth order Runge-Kutta algorithm with integrating time step dt = 10 -6 .

Contrast enhancement and image inversion

The image to process with the nonlinear oscillator network is the weak contrasted image of figure 9.(a). Its histogram is restricted to the range [0; 0.2], which means that the maximum gray level of the image (0.2) is the initial condition of at least one oscillator of the network, while the minimum gray level of the image (0) is also the initial condition of at least one oscillator. Therefore, the pixels with initial gray level 0 and 0.2, oscillate with the phase difference δ(t) predicted by figure 7.(b). In particular, as explained in section II.B.2, their phase difference δ(t) can be null for the processing times t = 3.96 × 10 -4 , 1.81 × 10 -4 , 3.5 × 10 -4 , 5.21 × 10 -4 ; minimum for t = 1.4 × 10 -4 , 4.64 × 10 -4 , 1.47 × 10 -3 and maximum for t = 3 × 10 -4 , 6.29 × 10 -3 , 1.64 × 10 -3 . As shown in figure 9.(b), (d), (f ) and (h), the image goes through local minima of contrast at the processing times corresponding to the zeros of δ(t) . Furthermore, the processing times providing the local minima of δ(t) realize an image inversion with a growing contrast enhancement (Fig. 9.(c), (g) and (j)). Indeed, since the minima of δ(t) are negative, for these processing times the minimum of the initial image becomes the maximum of the filtered image and vice-versa. Lastly, the local maxima of δ(t) allow to achieve local maxima of contrast for the corresponding processing times (Figs. 9.(e), (i), (k)). Note that the best enhancement of contrast is attained at the processing time t opt for which δ(t) is maximum. The histogram of each filtered image in fig. 9 also reveals the temporal dynamic of the network. Indeed, the width of the image histogram is periodically increased and decreased, which involves that the contrast of the corresponding filtered image is periodically enhanced or reduced.

Another, interesting feature of the realized contrast is determined by the plot of the network response at the processing time t opt [START_REF] Morfu | Image processing with a cellular nonlinear network[END_REF]. Indeed, this curve also represents the gray level of the pixels of the filtered image versus their initial gray level. Therefore, the horizontal axis corresponds to the initial gray scale, namely [0; 0.2], whereas the vertical axis represents

W i 0 W i (a) 0 0 . 1 0 . 2 0 2 . 5 5 W i 0 W i (b) o p t 0 0 . 1 0 . 2 0 2 . 5 5 W i 0 W i (c) 0 0 . 1 0 . 2 0 2 . 5 5 W i 0 W i (d) FIGURE 10
. Response of the nonlinear system for different nonlinearity parameters α at the corresponding optimal time topt (solide line) compared to a uniform rescaling (dotted line). The curves are obtained with eqs. ( 20) and ( 21) setting the time to the optimum value defined by the maximum of δ(t) (see fig. 8). In addition, we let the initial conditions W 0 i varying in the range [0; 0.2] in eqs. ( 20) and ( 21). (topt = 0.95 × 10 -3 ; α = 1.63). ω 0 = 1 the gray scale of the processed image. Such curves are plotted in figure 10 for different values of the nonlinearity parameter α, and at the optimal time defined by the maximum of δ(t). In fact, these times were established in section II.B.2 at figure 8. Moreover, to compare our nonlinear contrast enhancement to a uniform one, we have superimposed (in dotted line) the curve resulting from a simple multiplication of the initial gray scale by a scale factor. In fig 10.(a), since the response of the system for the lowest value of α is most often above the dotted lines, the filtered image at the processing time t opt = 1.75 × 10 -3 for α = 0.4 will be brighter than the image obtained with a simple rescaling.

As shown in Fig. 10. (b), increasing the nonlinearity parameter α to 1.05 involves an optimum time 1.66 × 10 -3 and symmetrically enhances the light and dark gray levels. When the nonlinearity parameter is adjusted to provide the greatest potential barrier (Figs. 10. (c) and (d)), the contrast of the medium gray level is unchanged compared to a simple rescaling. Moreover, the dark and light grays are strongly enhanced with a greater distortion when the potential barrier is maximum, that is for the greatest value of α (Fig. 10. (d)).

Gray level extraction

Considering processing times exceeding the optimal time t opt , we propose to perform a gray level extraction of the continuous gray scale represented figure 11.(a) [START_REF] Morfu | Image processing with a cellular nonlinear network[END_REF]. For a sake of clarity, it is convenient to redefine the white level by 0.2 whereas the black level remains 0.

For the 9 specific times presented in figure 11, the response of the system displays a minimum which is successively reached for each level of the initial gray scale. Therefore, time acting as a discriminating parameter, an appropriate threshold filtering allows to extract all pixels with a gray level in a given range. Indeed, in figure 11, the simplest case of a constant threshold V th = 0.25 provides 9 ranges of gray at 9 closely different processing times, which constitutes a gray level extraction. Moreover, owing to the response of the system, the width of the extracted gray level ranges reduces in the light gray. Indeed, the range extracted in the dark gray for the processing time t = 3.33 × 10 -3 (fig. 11.(c)) is approximatively twice greater than the range extracted in the light gray for t = 3.51 × 10 -3 (figure 11.(i)). To perform a perfect gray level extraction, the threshold has to match with a slight offset the temporal evolution of the minimum attained by the response of the system. Under these conditions, the width of the extracted gray range is set by the value of this offset.

Note that, the response of the system after the optimal processing times also allows to consecutively enhance fragment of the image with different levels of brightness, which is also an important feature of the image processing. For instance, in Belousov-Zhabotinsky-type media this property of the system enabled Rambidi and co-workers to restore individual components of the picture when the components are overlapped [START_REF] Rambidi | Image processing using light-sensitive chemical waves[END_REF]. Therefore, we trust that considering the temporal evolution of the image loaded in our network could give rise to other interesting image processing operations.

Image encryption

Cryptography is another field of application of nonlinear systems. In fact, the chaotic behavior of nonlinear systems can sometimes produce chaotic-like waveforms which can be used to encrypt signals for secure communications [START_REF] Cuomo | Circuit Implementation of Synchronized chaos with applications to communications[END_REF][START_REF] Dedieu | Chaos shift keying: modulation and demodulation of a chaotic carrier using self-synchronizing Chua's circuits[END_REF]. Even if many attempts to break the encryption key of these cryptosystems and to retrieve the information have been reported [START_REF] Short | Unmasking a hyperchaotic communication scheme[END_REF][START_REF] Udaltsov | Cracking chaos-based encryption systems ruled by nonlinear time delay differential equations[END_REF], cryptography based on the properties of chaotic oscillators still attracts the attention of researchers owing to the promising applications of chaos in data transmission field [START_REF] Kwok | A fast image encryption system based on chaotic maps with finite precision representation[END_REF]. Contrary to most studies, where the dynamics of a single element is usually considered, we propose here a strategy of encryption based on the dynamics of a chain of nonlinear oscillators. More precisely, we consider the case of a noisy image loaded as initial condition in the inertia network introduced in section II.B. In addition, we add a uniform noise over [-0.1; 0.1] to the weak contrasted picture of the coliseum represented in Fig. 9.(a). Since the pixels of the noisy image take a gray level in the range [-0.1; 0.3], an appropriate change of scale is realized to reset the dynamics of the gray levels to [0; 0.2]. The resulting image is then loaded as initial condition in the network. For a sake of clarity, the filtered images are presented at different processing times with the corresponding system response in figure 12.

0 0 . 2 (a) 0 0 . 2 0 6 W i 0 W i (b) 0 0 . 2 0 6 W i 0 W i (c) 0 0 . 2 0 6 W i W i (d) 0 0 . 2 0 6 W i 0 W i (e) 0 0 . 2 0 6 W i 0 W i (f) 0 0 . 2 0 6 W i W i (g) 0 0 . 2 0 6 W i 0 W i (h) 0 0 . 2 0 6 W i 0 W i (i) 0 0 . 2 0 6 W i W i (j)
Before the optimal time, we observe the behavior described in section III.A.1: the image goes through local minima and maxima of contrast until the optimum time t opt = 1.64 × 10 -3 , where the best contrast enhancement is realized (Fig 12 .(a)).

Next, for processing times exceeding t opt , the noisy part of the image seems to be amplified while the coherent part of the image begins to be less and less perceptible (see Fig. 12.(b) and 12.(c) obtained for t = 3.28 × 10 -3 and t = 6.56 × 10 -3 ). Lastly, for greater processing times, namely t = 8.24 × 10 -3 and t = 9.84 × 10 -3 , the noise background has completely hidden the coliseum which constitutes an image encryption.

Note that this behavior can be explained with the response of the system, as represented below each filtered image of Fig. 12. Indeed, until the response of the system versus the initial condition does not display a "periodic-like" behavior, the coherent part of the image remains perceptible (Fig. 12.(a) and (b)). By contrast, as soon as a "periodicity" appears in the system response, the coherent image begins to disappear (Fig. 12.(c)). Indeed, the response of Fig. 12.(c) shows that 4 pixels of the initial image with 4 different gray levels take the same final value in the encrypted image (see the arrow lines). Therefore, the details of the initial image, which corresponds to quasi-uniform area of the coherent image, are merged and thus disappear in the encrypted image. Despite the previous merging of gray levels, since noise induces sudden changes in the gray levels of the initial image, the noise conserves its random feature in the encrypted image. Moreover, since the system tends to enlarge the range of amplitude, the weak initial amount of noise is strongly amplified whatever the processing time exceeding t opt . The periodicity of the system response can then be increased for larger processing times until only the noisy part of the image is perceptible (Fig. 12.(d) and (e)). A perfect image encryption is then realized.

To take advantage of this phenomenon for image encryption, the coherent information, that is the enhanced image of Fig. 12.(a), must be restored using the encrypted image of Fig. 12.(e). Fortunately, owing to the absence of dissipation, the nonlinear systems is conservative and reversible. It is thus possible to go back to the optimal time that is when the information was the most perceptible.

However, the knowledge of the encrypted image is not sufficient to completely restore the coherent information, since at the encryption time, the velocity of the oscillators was not null. Consequently, it is necessary to know both the position and the velocity of all particles of the network at the time

W i m -a A D 6 3 3 J N -K ( W i -m ) ( W i -m + a ) / 1 0 m + a A D 6 3 3 J N -K P ( W i ) = ( W i -m + a ) ( W i -m ) ( W i -m -a ) K 2 / 1 0 0 W i = P ( W i ) m V T + W i 0 FIGURE 13
. Sketch of the elementary cell of the inertial system. m and α are adjusted with external DC sources whereas -K is the inverting amplifier gain obtained using T L081CN operational amplifier. The 1N 4148 diode allows to introduce the initial conditions W 0 i .

of encryption. The information can then be restored solving numerically eq. ( 23) with a negative integrating time step dt = -10 -6 . Under these conditions, the time of encryption constitutes the encryption key.

B. Electronic implementation

The elementary cell of the purely inertial systems can be developed according to the principle of figure 13 [START_REF] Morfu | Experiments on an electrical nonlinear oscillators network[END_REF]. First, a polynomial source is realized with analog AD633JN Z multipliers and classical inverting amplifier with gain -K. Taking into account the scale factor 1/10 V -1 of the multipliers, the response of the nonlinear circuit to an input voltage W i is given by

P (W i ) = K 2 100 (W i -m)(W i -m -α)(W i -m + α), (24)
where the roots m, m -α, m + α of the polynomial circuit are set with three different external DC sources. As shown in figure 14, the experimental characteristic of the nonlinear source is then in perfect agreement with its theoretical cubic law (24).

1 . 5 2 2 . 5 3 3 . 5 4 -0 . 4 -0 . 2 0 0 . 2 0 . 4 W i ( V o l t ) P ( W i ) ( V o l t )
Next, a feedback between the input/output of the nonlinear circuits is ensured by a double integrator with time constant RC such that

W = - K 2 100R 2 C 2 (W i -m + α)(W i -m -α)(W i -m)dt. ( 25
)
Deriving twice eq. ( 25), the voltage W i at the input of the nonlinear circuit obeys to

d 2 W i dt 2 = - K 2 100R 2 C 2 (W i -m + α)(W i -m -α)(W i -m), (26)
which corresponds exactly to the equation of the purely inertial system (13) with

ω 0 = K/(10RC). (27)
Lastly, the initial condition W 0 i is applied to the elementary cell via a 1N 4148 diode with threshold voltage V T = 0.7 V . Indeed, we adjust the diode anode potential to W 0 i + V T with an external DC source involving that the diode cathode potential is initially set to W 0 i . Then, according to section III, the circuit begins to oscillate in the range [W 0 i ; 2m -W 0 i ], while the potential of the diode anode remains V T + W 0 i . Assuming that m > W 0 i /2, which is the case of our experiments, the diode is instantaneously blocked once the initial condition is introduced. Note that using a diode to set the initial condition presents the main advantage to "balance" the effect of dissipation inherent in electronic devices. Indeed, the intrinsic dissipation of the experiments tends to reduce the amplitude of the oscillations W 0 i . As soon as the potential of the diode cathode is below W 0 i , the diode conducts instantaneously, introducing periodically the same initial condition in the elementary cell. Therefore, the switch between the two states of the diode presents the advantage to refresh the oscillations amplitude to their natural worths as in absence of dissipation.

In summary, the oscillations are available at the diode cathode and are represented in figure 15

W i ( V o l t ) W i ( t o p t ) ( V o l t ) 0 FIGURE 16
. Response of the system to a set of initial conditions W 0 i ∈ [0; 0.2] at the optimal time. The solid line is obtained with eqs. ( 20), ( 21) and ( 27 1V is obtained corresponding to the phase opposition W 1 (t opt ) = 0V and W 2 (t opt ) = 5.1V . Therefore, the weak difference of initial conditions between the oscillators is strongly increased at the optimal time t opt . Despite a slight discrepancy of 11% for the optimal time, mainly imputable to the component uncertainties, a purely inertial nonlinear system is then implemented with the properties of section III.

To perfectly characterize the experimental device, we now focus on the response of the nonlinear system to different initial conditions in the range [0V ; 0.2V ]. The plot of the voltage reached at the optimal time t opt = 1.46ms versus the initial condition is compared in figure 16 to the theoretical curve obtained for the optimum time defined in section II.B.2, namely 1.64 ms. The experimental response of the system is then qualitatively confirmed by the theoretical predictions, which allows to valid the experimental elementary cell for the contrast enhancement presented in section III.A.1.

Lastly, we also propose to investigate the response of the system after the optimum time, since it allows the extraction of gray levels. In order to enhance the measures accuracy, we extend the range of initial conditions to [0, 0.5V ] instead of [0, 0.2V ]. The corresponding experimental optimal time becomes t opt = 564µs, whereas the theoretical ones, deduced with the We have also reported the experimental device response for three different times beyond the optimal time t opt = 564µs in figure 17.(b), (c), (d), namely for the experimental times t = 610µs, t = 675µs and t = 720µs. Since a time scale factor 610/564 = 1.1684 exists between the experimental and the theoretical optimal time, we apply this scale factor to the three previous experimental times. It provides the theoretical times 713µs, 789µs, 841µs.

W i ( V o l t ) W i ( t o p t ) ( V o l t ) -0 . 1 0 0 . 1 0 . 2 0 . 3 0 . 4 0 . 5 0 . 6 -1 0 1 2 3 4 5 6 0 (a) -0 . 1 0 0 . 1 0 . 2 0 . 3 0 . 4 0 . 5 0 . 6 -0 . 5 0 0 . 5 1 1 . 5 2 2 . 5 3 3 . 5 4 4 . 5 W i ( V o l t ) W i ( t o p t ) ( V o l t ) 0 (b) 0 0 . 5 1 1 . 5 2 2 . 5 3 3 . 5 4 4 . 5 5 W i ( V o l t ) -0 . 1 0 0 . 1 0 . 2 0 . 3 0 . 4 0 . 5 0 . 6 W i ( t o p t ) ( V o l t ) 0 (c) W i ( t o p t ) ( V o l t )
For each of these 3 times, we are then able to compare the experimental response to the theoretical one deduced by letting the initial condition vary in [0; 0.5V ] in eqs. ( 20), ( 21) and ( 27). Despite some slight discrepancies, the behaviour of the experimental device is in good agreement with the theoretical response of the system for the three processing times exceeding the optimal time. Therefore, the extraction of gray levels, presented in section III.A.2, is electronically implemented with this elementary cell.

IV. REACTION-DIFFUSION SYSTEMS

A. one dimensional lattice

The motion equation (4) of the nonlinear mechanical chain can also describe the evolution of the voltage at the nodes of a nonlinear electrical lattice. This section is devoted to the presentation of this nonlinear electrical lattice.

The nonlinear lattice is realized by coupling elementary cells with linear resistors R according to the principle of figure 18.(a). Each elementary cell

R R R R U n U n -1 U n + 1 C R N L R N L C R N L C (
= 3.078KΩ, V b = 1.12V , Va = 0.545V , β = 1.
consists of a linear capacitor C in parallel with a nonlinear resistor whose current-voltage characteristic obeys to the cubic law

I N L (u) = βu(u -V a )(u -V b )/(R 0 V a V b ), (28)
where 0 < V a < V b are two voltages, β is a constant and R 0 is analog to a weighting resistor.

The nonlinear resistor can be developed according to two different methodologies. The first way to obtain a cubic current is to consider the circuit of figure 18.(b) with three branches [START_REF] Comte | Etude d'une ligne non linéaire de type Nagumo-Neumann[END_REF][START_REF] Binczak | Experiments on a nonlinear electrical reaction-diffusion line[END_REF]. A linear resistor R 3 , a negative resistor and another linear resistor R 1 are successively added in parallel thanks to 1N 4148 diodes. Due to the switch of the diodes, the experimental current-voltage characteristic of figure 19 asymptotically displays a piecewise linear behaviour with successively a positive slope, a negative one and lastly a positive one.

This piecewise linear characteristics is compared to the cubic law (28) which presents the same roots V a , V b and 0 but also the same area below the characteristic between 0 and V a . This last conditions leads to β = 1 and R 0 = 3.078 KΩ (Morfu, 2002c). An alternative way to realize a perfect cubic nonlinear current is to use P o l y n o m i a l g e n e r a t i o n c i r c u i t s (Comte and Marquié, 2003).

R 0 P ( U ) U R N L I N L ( U ) U I N L ( U )
(u) = βu(u -V a )(u -V b )/(V a V b ) + u as shown in figure 20
This polynomial voltage is realized with AD633JN Z multipliers and classical T L081CN operational amplifiers. A resistor R 0 ensures a feedback between the input/output of the nonlinear source such that the Ohm's law applied to R 0 corresponds to the cubic current (28):

P (u) -u R 0 = I N L (u). (29)
As shown in Fig. 21, this second method gives a better agreement with the theoretical cubic law (28).

Applying the Kirchhoff laws, the voltage U n at the n th node of the lattice obeys to involves for n = 1 and n = N respectively

C dU n dτ = 1 R (U n+1 + U n-1 -2U n ) -I N L (U n ), (30
C dU 1 dτ = 1 R (U 2 -U 1 ) -I N L (U 1 ), (31)
C dU N dτ = 1 R (U N -1 -U N ) -I N L (U N ). (32)
Next, introducing the transformations

W n = U n V b , D = R 0 R αβ, t = τ R 0 αCβ , (33)
yields the discrete Nagumo equation in its normalized form,

dW n dt = D(W n+1 + W n-1 -2W n ) + f (W n ). (34)
Therefore, an electronic implementation of the overdamped network presented in section II.A is realized.

B. Noise filtering of a 1D signal

One of the most important problem in signal or image processing is probably to remove noise from a coherent information. In this section, we intensively develop the principle of nonlinear noise filtering inspired by the overdamped systems [START_REF] Marquié | Diffusion effects in a nonlinear electrical lattice[END_REF]. In addition, using the electrical nonlinear network introduced in section IV.A, we also present an electronic implementation of the filtering tasks .

Theoretical analysis

To investigate the response of the overdamped network to a noisy signal loaded as initial condition, we first consider the simple case of a constant signal with a sudden change of amplitude. Therefore, we study the discrete normalized Nagumo equation

dW n dt = D(W n+1 + W n-1 -2W n ) + f (W n ), (35)
with

f (W n ) = -W n (W n -α)(W n -1) in the specific case α = 1/2.
Furthermore, the initial condition applied to the cell n is assumed to be uniform for all cells, except for the cell N/2 where a constant perturbation b 0 is added; namely:

W n (t = 0) = V 0 ∀n = N 2 W N/2 (t = 0) = V 0 + b 0 . ( 36
)
The solution of eq. ( 35) to the initial condition (36) can be expressed under the following form

W n (t) = V n (t) + ǫb n (t) (37)
Inserting eq. ( 37) in eq. ( 35), we collect the terms of order 0 and 1 in ǫ with the reductive perturbation methods to obtain the set of differential equations [START_REF] Taniuti | Reductive perturbation method in nonlinear wave propagation[END_REF][START_REF] Taniuti | Perturbation Method for a Nonlinear Wave Modulation[END_REF]:

dV n dt = D(V n+1 + V n-1 -2V n ) + f (V n ) (38) db n dt = D(b n+1 + b n-1 -2b n ) -(3V 2 n -2V n (1 + α) + α)b n (39)
Assuming that V n is a slow variable, eq. ( 38) reduces to

dV n dt = f (V n ), (40)
which provides the response of the system to a uniform initial condition V 0 (see details in Appendix A):

V (t) = 1 2 1 + V 0 -1 2 (V 0 -1 2 ) 2 -V 0 (V 0 -1)e -t 2 . ( 41
)
Next, to determine the evolution of the additive perturbation, it is convenient to consider a perturbation under the following form

b n (t) = I n (2Dt)g(t), (42)
where I n is the modified Bessel functions of order n [START_REF] Abramowitz | Handbook of Mathematical Functions[END_REF]. Substituting eq. ( 42) in eq. ( 39), and using the property of the modified Bessel function:

dI n (2Dt) dt = D(I n+1 + I n-1 ), (43)
we obtain straightforwardly

dg dt = -2Dg -3V 2 n -2V n (1 -α) + α g, (44)
that is

dg g = -2Ddt -3V 2 n -2V n (1 -α) + α dt. (45)
Noting that

df (V n ) dt = -3V 2 n -2V n (1 -α) + α dV n dt , (46)
and deriving eq. ( 40) versus time, we get

V ′′ n V ′ n = -3V 2 n -2V n (1 -α) + α , (47)
where V ′ n et V ′′ n denote the first and second derivative versus time. Combining eq. ( 47) and eq. ( 45) allows to express g(t) as:

g(t) = Ke -2Dt dV n dt ( 48
)
where K is an integrating constant. Deriving eq. ( 41), we obtain g(t) and thus the evolution of the perturbation:

b n (t) = K I n (2Dt)e -2Dt e -t/2 8 V 0 (V 0 -1 2 )(V 0 -1) (V 0 -1 2 ) 2 -V 0 (V 0 -1)e -t/2 3/2 (49)
Writing b n (t = 0) = b 0 n , provides the value of the integrating constant K. The evolution of the perturbation b n (t) is then ruled by:

b n (t) = b 0 n 8 I n (2Dt)e -2Dt e -t 2 (V 0 -1 2 ) 2 -V 0 (V 0 -1)e -t 2 3 2 . ( 50
)
Lastly, in the case of multiple perturbations, the perturbation at the n th node of the lattice obeys to

b n (t) = n ′ b 0 n ′ 8 I n ′ -n (2Dt)e -2Dt e -t 2 (V 0 -1 2 ) 2 -V 0 (V 0 -1)e -t 2 3 2 , ( 51
)
where I n ′ -n is the modified Bessel function of order n ′ -n.

Eq. ( 41) shows that the evolution of the constant background does not depends on the coupling D. By contrast, eq. ( 51) reveals that the coupling D can be tuned to speed-up the diffusion of the perturbation without affecting the constant background. Therefore, in signal processing context, this property can be used to develop a noise filtering tool.

Theoretical and numerical results

In order to valid the theoretical analysis developed in section IV.B.1, we have solved numerically eq. ( 35) using a fourth order Runge-Kutta algorithm with an integrating time step dt = 10 -3 . Moreover, a uniform initial condition V 0 = 0.4 is loaded for all the N = 48 cells of the network excepted for 41) and ( 51), (o) signs: numerical results.

the 24 th cell. Indeed, for this cell, an additive perturbation b 0 = 0.2 is superimposed onto the constant background V 0 in order to match exactly the initial condition (36) considered in the theoretical section IV.B.1.

We have first investigated the evolution of both the constant background and the perturbation versus time. In fig. 22, the numerical results plotted with (o) signs match with a perfect agreement the theoretical results predicted by eqs. ( 41) and (51). Moreover, the curves (a) of Fig. (22) shows that the constant background given by eq. ( 41) is unaffected by the nonlinear systems whatever the coupling value D. By contrast, the behaviour of the system for the additive perturbation b 0 depends on the coupling parameter D (curves (b) and (c)). Indeed, for weak coupling values, namely D = 0.5, the perturbation slowly decreases and seems to be quasi-unchanged, whereas for D = 5, the curve (c) exhibits a greater decreasing behavior. After the time t = 0.4, the perturbation is significantly reduced for D = 5. Therefore, the coupling parameter D can be tuned to speed up the diffusion of the perturbation without disturbing the constant background. Furthermore, the time acts as a parameter which adjusts the filtering of the perturbation. The state of the lattice for two different processing times is represented in fig. 23 ( shows that the perturbation is not filtered for D = 0.5 and for the same processing time t = 2. Furthermore, in both cases the constant background is slowly attracted by the nearest stable state, that is 0 in our case. Note that the spatiotemporal views of fig. 24 also reveal that the noise filtering is performed for D = 5 and a processing time t = 2. Lastly, to valid the processing task realized by the overdamped system, we propose to remove the noise from a more complex signal, that is a noisy sinusoïdal signal. The signal is first sampled with a total number of samples corresponding to the size of the overdamped network, namely N . Next, a serial to parallel conversion is realized to load the N samples at the nodes of the 1D lattice. Therefore, we are led to consider the distribution of initial

x n = A cos 2π 2n N + 1 2 + η n , ( 52
)
where η n is a discrete white gaussian noise of RM S amplitude σ = 0.15.

A and 2/N represent respectively the amplitude and the frequency of the coherent signal. We first numerically investigate the response of the network with the coupling D = 0.5. As in the case of a constant background corrupted by a local perturbation, the system is unable to remove the noise from the sinusoidal signal for both the processing times presented in fig. 25.(b) and (d). By contrast, for the favorable value of the coupling D = 5, the noise is completely filtered at the processing time t = 1 as shown in fig. 25.(e).

Experimental results

To valid the electronic implementation of the nonlinear noise filtering tool, we consider the nonlinear electrical lattice introduced in section IV.A with the nonlinear resistor of figure 18.(b). In order to match the coupling value D = 5 and D = 0.5, the coupling resistor R is set to R = 300Ω and R = 3KΩ respectively. Moreover, all results are presented in normalized units using the transformation (33) to allow a direct comparison with the theoretical analysis of section IV.B.2. First, we experimentally report in fig 26

the temporal evolution of the set of initial conditions consisting of a constant signal locally corrupted by a perturbation. As predicted in the theoretical section, the constant background is unaffected whatever the coupling value (curve (a)) whereas when the coupling is adjusted to its favorable value D = 5, the perturbation can be removed after a normalized processing time t = 0.4 (curve (c)). This results is also confirmed by the spatial response of the system at two different processing times. Indeed, as shown in Fig. 27, the state of the lattice for t = 2 and t = 4 provides the signal without the perturbation only if the coupling D is chosen equal to 5. Lastly, we propose to filter the noisy sinusoidal signal of fig. 28.(a). After a processing time t = 0.6, the noise is completely removed for the coupling D = 5 (fig. 28.(c)), which is not the case if the coupling is set to D = 0.5 (fig. 28.(b)). Therefore, with a suitable choice of both processing time and resistor coupling, a noise filtering tool inspired by the properties of the nonlinear overdamped network is electronically implemented. Moreover, according to the transformation (33), the processing time could be adjusted by the worth of the capacitor C to match real time processing constraints.

C. 2D filtering: Image processing

We now numerically extend the properties of the one dimensional lattice to a two dimensional network. We are led to consider a CN N whose cell state W i,j , representing the gray level of the pixel number i, j, obeys to the following set of equations:

dW i,j dt = f (W i,j ) + D (k,l)∈N r (W k,l -W i,j ), i, j = 2...N -1, 2...M -1, ( 53
)
where N r = {(i -1; j), (i + 1, j), (i, j + 1), (i, j -1)} is the set of the four nearest neighbors, N × M the image size and f (W i,j ) represents the nonlinearity. The boundary conditions for the edges of the image express

dW 1,j dt = f (W 1,j ) + D(W 1,j-1 + W 2,j + W 1,j+1 -3W 1,j ), j = 2..M -1 dW N,j dt = f (W N,j ) + D(W N,j-1 + W N -1,j + W N,j+1 -3W N,j ), j = 2..M -1 dW i,1 dt = f (W i,1 ) + D(W i-1,1 + W i+1,1 + W i,2 -3W i,1 ), i = 2..N -1 dW i,M dt = f (W i,M ) + D(W i-1,M + W i+1,M + W i,M-1 -3W i,M ), i = 2..N -1
while for the image corners, we consider the two nearest neighbors, that is

dW 1,1 dt = f (W 1,1 ) + D(W 2,1 + W 1,2 -2W 1,1 ), dW N,M dt = f (W N,M ) + D(W N,M-1 + W N -1,M -2W N,M ), dW N,1 dt = f (W N,1 ) + D(W N -1,1 + W N,2 -2W N,1 ), dW 1,M dt = f (W 1,M ) + D(W 2,M + W 1,M-1 -2W 1,M ).

Noise filtering

The initial condition applied to the cell i, j of the network corresponds to the initial gray level W 0 i,j of the noisy image represented in fig. 29. The image after a processing time t is obtained noting the state W i,j (t) of all cells of the network at this specific time t [START_REF] Comte | Noise removal using a twodimensional diffusion network[END_REF].

In fig. 30, we have reported the filtered image obtained at the processing times t = 1, t = 3, t = 6, t = 9 and for the coupling values D = 0.075, D = 0.1, D = 0.2 and D = 0.3 respectively. The bistable behaviour of the system established in section II.A.1 involves a natural evolution of the image towards the two stable states of the system, namely 0 and 1. It means that, as time grows, the image evolves to a black and white pattern. Therefore, to realize a correct noise filtering, the coupling parameter and the processing time must be adjusted.

For the lowest coupling value D = 0.075, fig. 30 shows that the noise is not removed before the image is binarized. For the coupling parameter D = 0.2 and D = 0.3, even if the noise is quickly removed, the filtered image becomes blurred for t = 6 and t = 9 (Fig. 30.(k),(l), (o),(p)). Therefore, these settings of the coupling parameter are unappropriate. In fact, in fig 30 .(f ) and (g), we have obtained the filtered image with the best setting of the coupling and of the processing time, that is a coupling D = 0.1 and the processing times t = 3 or t = 6. Indeed, the filtered images are neither blurred nor binarized. Moreover, the system does not only remove the noise, it also enhances the contrast of the initial image.

Edge Filtering

Because of a strong relationship between edge and object recognition, edge detection constitutes one of the most important steps for image recognition. Indeed, scene information can often be interpreted thanks to the edges. Classical edge detection algorithms are based on a second order local derivative operator [START_REF] Gonzalez | Digital image processing[END_REF] while nonlinear techniques of edge enhancement are mainly inspired by the properties of reaction-diffusion media [START_REF] Rambidi | Image processing using light-sensitive chemical waves[END_REF]Chua and Yang, 1988).

We propose here a strategy of edge detection based on the propagation properties of the nonlinear diffusive medium (Comte et al, 2001). The image loaded in the 2-Dimensional network is the black and white picture of fig.

31.(a).

We have established in section II.A.2 that a 1D lattice modeled by the Nagumo equation supports kink and anti-kink propagation owing to the bistable nature of the nonlinearity. Indeed, if the nonlinearity threshold parameter α < 1/2, the stable state 1 propagates, while if α > 1/2 the stable state 0 propagates. Therefore, extending this property to a 2 Dimensional network allows to perform either erosion for α > 1/2 or dilation for α < 1/2, which are basic mathematical morphology operations, commonly performed in image processing [START_REF] Serra | Introduction to mathematical morphology[END_REF]. Moreover, if the initial image is subtracted to the image obtained with the network obeying to eq. ( 53) , we can deduce after a processing time t, the contours of the image. We have reported in Fig. 31.(b) the contour of a black and white image and its profile obtained with this method. The profile of the contour shows that its resolution is about 10 pixels, which is too important to allow a good edge enhancement of a more complex image. This poor resolution is mainly imputable to the spatial expansion of the kink which results from the initial condition loaded in the lattice. Since the kink expansion reduces with the coupling, a natural solution consists in lowering the coupling. Unfortunately, the existence of the propagation failure effect provides a lower bound of the coupling D * and thus hinder a contour detection with a good resolution. An alternative solution can be developed by using a nonlinearity which forbids the existence of the propagation failure effect. Indeed, it has been shown for dissipative media [START_REF] Bressloff | Exact travelling wave solutions of an "integrable" discrete reaction-diffusion equation[END_REF] or for systems where both inertia and dissipation are taken into account [START_REF] Comte | Dissipative lattice model with exact travelling discrete kink-soliton solutions: Discrete breather generation and reaction diffusion regime[END_REF], that an inverse method allows to define a nonlinear function for which exact discrete propagative kinks exists. Especially, in the purely dissipative case, such function expresses where ǫ = 0.5, a 2 = 0.9, a 0 = 1.483 and a 1 = -0.742 to ensure that the zeros of the nonlinearity remain 0, 1/3 and 1. As expected, when this new nonlinearity is numerically implemented, the resolution of the detected contour in Fig 31 .(c) reduces to 3 pixels. Note that edge enhancement with the nonlinear overdamped network is not restricted to a black and white image. Indeed, the concept is based on the propagation properties of the system and can be extended to the case of an image with 256 gray levels. For instance, we numerically propose to show the contour enhancement of Fig. 32.(a) by considering the methodology used for the edge detection of the black and white picture. The simulation results are summarized in Fig. 32 for different processing times in the favorable case of the nonlinear function (54). It is clear that once again the time allows to adjust the quality of the processing. Indeed, for processing times below t = 1, the edges of the image details are not revealed, whereas for processing times exceeding 1.33, the details begin to disappear. Furthermore, as times grows, the contours of the image are less and less thin owing to the propagation mechanism. The best contour enhancement is thus performed when the image details have not yet disappeared and when the enhanced contours remain sufficiently thin. In fact, this situation corresponds to the intermediate processing time t = 1.33 (Fig. 32.(e)).

(a) (b) (c) (d) (e) (f) (g) (h) (i) (j) (k) (l) (m) (n) (o) (p)
f (W i,j ) = Dǫ (1 -a 2 /2) -(a 0 W i,j + a 1 ) 2 - Da 2 (a 0 W i,j + a 1 ) 1 -(a 0 W i,j + a 1 ) 2 + 2D(a 0 W i,j + a 1 ), (54)

Extraction of regions of interest

As explained in the previous subsections, in the case of the cubic nonlinearity, a nonlinearity threshold α = 0.5 allows to perform noise filtering, while considering α = 0.5 provides the contour of an image with a poor resolution. Moreover, the nonlinearity f (W ) can be determined using an inverse method to optimize the filtering task. Therefore, the choice of the nonlinearity is of crucial interest to develop interesting and powerful image processing tools. In this section, we go one step further by proposing a new nonlinearity to extract the regions of interest of an image representing the soldering between two rods of metal [START_REF] Morfu | Experiments on an electrical nonlinear oscillators network[END_REF].

The noisy and weakly contrasted image of Fig. 33 presents 4 regions of interest:

• First, the two rods of metal constitute the background of the image in light gray; • The stripe in medium gray at the center of the image represents the "soldered joint"; • A white spot corresponds to a "projection" of metal occurring during the soldering of the two rods of metal; • A dark gray spot represents a gaseous inclusion inside the soldering joint.

(a) (b) (c) (d) (e) (f) (g) (h) (i) (j) FIGURE 32.
Contour enhancement of a an image with 256 gray levels realized with the modified nonlinearity (54). (a) Initial image. (b), (c), (d), (e), (f ), (g), (h), (i) and (j) Filtered image for the respective processing times t = 0.33, t = 0.66, t = 1, t = 1.33, t = 1.66, t = 2, t = 2.33, t = 2.66 and t = 3.

Limit of the bistable Network

We first discuss the inability of the bistable overdamped network ruled by eq. ( 53) to extract the 4 objects of the image. As explained in section II.A.1, the bistability is ensured by using the cubic nonlinearity (3). According to the mechanical description of the bistable system presented in section II, a pixel of the image is analog to a particle experiencing a double well potential φ(W ) = -

D (k,l)∈N r W k,l -W i,j , (55)
and the nonlinear force f (W i,j ). This property of the system involves that for sufficiently large time, the network is organized near the two stable states set by the nonlinearity, namely 0 and 1. In image processing context, it means that the resulting filtered image tends to an almost black and white pattern. Fig. 35 confirms this evolution of the filtered image versus the processing time since, when a cubic nonlinearity is considered, a quasi-black and white image is obtained at the time t = 3000 (Fig. 35.(c)).

Note that for none of the proposed processing times, the bistable system was able to properly remove the noise and to enhance the contrast of the regions of interest. Indeed, for t = 4 the noise is reduced but the details of the image begin to disappear (Fig. 35.(a)). Especially, the projection is merged into the background for t = 10, revealing that the bistable nature of the system destroys the coherent information of the initial image (Fig. 35.(b)). Therefore, the inability of the overdamped system to extract the regions of interest is directly imputable to the bistable nonlinear force f (W ).

The multistable network

To solve this problem and to keep the coherent structure of the image, we introduce a nonlinearity with a multistable behavior. For instance, the following nonlinear force The multistable behavior of the network obeying to eq. ( 53) with the sinusoidal force (56) can be established by considering the uncoupled case. Setting D = 0 in eq. ( 53), we obtain

f (W ) = -β(n -1) sin 2π(n -1)W , ( 56
) W i , j F ( W i , j ) 0 0 . 5 1 -0 . 0 2 0 0 . 0 2 b / p
dW i,j dt = -β(n -1) sin 2π(n -1)W i,j . (57)
The stability analysis of the system can be performed with the methodology developed in section II.A.1 by considering the roots of the sinusoidal force (56). According to the sign of the derivative of the sinusoidal force, we can straightforwardly deduce that the unstable steady states of the system are given by

W thk = (2k + 1)/(2(n -1)) with k ∈ Z, (58)
while the stable steady states are defined by

W * k = k/(n -1) with k ∈ Z. (59)
The eq. ( 57) is solved in Appendix C to give the temporal evolution of an overdamped particle experiencing the multistable potential of fig. 36 in the uncoupled case. If k denotes the nearest integer of (n -1)W 0 i,j , and W 0 i,j the initial position of the particle, the displacement W i,j (t) of the particle obeys to

W i,j (t) = 1 π(n -1) arctan tan(π(n -1)W 0 i,j )e -β(n-1) 2 2πt + k n -1 . (60)
The multistable behaviour of the system is illustrated in fig. 37, where we have reported the temporal evolution of a particle submitted to different initial conditions in the range [0; 1]. It is clear that the unstable steady states of the system W thk act as thresholds, while the stable steady states W * k correspond to attractors. Indeed, the final state of the particle depends on the value of the initial condition compared to the thresholds W thk . Especially, if we neglect the transient, the asymptotic behavior of the uncoupled network reduces to the following rules

if 2k -1 2(n -1) < W 0 i,j < 2k + 1 2(n -1) W i,j (t → +∞) = k (n -1)
.

Therefore, the asymptotic functioning (61) of the uncoupled network proves the multistable behavior of the system. We now numerically use this multistable feature to extract the regions of interest of the image. In the coupled case, a pixel with initial gray level W 0 i,j can take one of the n possible stable states according to the competition between the sinusoidal force and the resulting elastic force. The specific case n = 5 is numerically reported in Fig. 38. Unlike the bistable network, the noise is quickly removed without disturbing the coherent structure of the image consisting of "the projection", "the gaseous inclusion", "the background" and the "soldered joint" (Fig. 38 ). An extraction of the interest regions of the image is then performed with this overdamped network.

Electronic implementation of the multistable network

The electronic implementation of the multistable network is realized according to the methodology of Fig. 39, by coupling elementary cells with linear resistors. Each elementary cell includes a capacitor in parallel with a nonlinear resistor

j -1 j j + 1 i -1 i i + 1 R N L C R R R R R R R I N L ( U i , j ) U i , j U i -1 , j + 1 U i -1 , j U i -1 , j -1 U i + 1 , j + 1 U i + 1 , j U i + 1 , j -1 U i , j + 1 U i , j -1 R R R FIGURE 39
. Electronic sketch of the multistable nonlinear network. R Represents the coupling resistor, C a capacitor and R NL a nonlinear resistor.

whose current-voltage characteristics can be approximated by the sinusoïdal law on the range [-2V ; 2V ]:

I N L (U ) ≃ I M sin(2πU ). (62)
The methodology detailed in the section IV.A to realize the cubic nonlinearity with a polynomial source can be once again used to obtain the sinusoidal law (62). First, a least square method at the order 15 allows to fit the sinusoidal expression (62) by a polynomial law P (U ) in the range [-2V ; 2V ]. It provides the coefficients of the polynomial source P (U ) to generate the sinusoidal current

I N L (U ) = P (U )/R 0 . ( 63
)
The experimental current-voltage characteristics is compared in of the least square method. However, enhancing the agreement with the sinusoidal law presents the main disadvantage to considerably increase the number of electronic components used for the realization of the nonlinear resistor. Nevertheless, at the order 15, the experimental nonlinear current presents 9 zeros and its derivative ensures the existence of 5 stable steady states and 4 unstable steady states. It is thus not of crucial interest to increase the order of the approximation, provided that the nonlinear resistor exhibits the multistability.

U ( V o l t ) U t h 1 U t h 2 U t h 3

Applying the Kirchhoff laws to the electrical network of Fig. 39, we deduce the differential equation, which rules the evolution of the voltage U i,j at the nodes (i, j)

C dU i,j dτ = -I N L (U i,j ) + 1 R (k,l)∈N r (U k,l -U i,j ). (64)
In eq. ( 65), N r = {(i; j -1), (i; j + 1), (i -1; j), (i + 1; j)} denotes the neighborhood and τ represents the experimental time.

Next, the transformations

τ = tR 0 C, β = I M R 0 (n -1) 2 , U i,j = W i,j (n -1) -2 and D = R 0 R , (65)
lead to the normalized equation

dW i,j dt = P (W i,j (n -1) -2) n -1 + D (k,l)∈N r (W k,l -W i,j ). ( 66
)
The normalization is completed by noting that for W i,j ∈ [0; 1], that is for U i,j ∈ [-2V ; 2V ],

P W i,j (n -1) -2 = -R 0 I N L W i,j (n -1) -2 ≃ -β(n -1) 2 sin(2π(n -1)W i,j ). ( 67
)
The experimental network described by eq. ( 66) appears as an analog simulation of the normalized multistable network used for image processing.

Let us finally reveal the multistable behaviour of the elementary cell of the experimental network by investigating its response to different initial conditions in the uncoupled case. In addition, to allow a direct comparison with the theoretical expression (60), all the results are presented in normalized units in Fig. 40.(a). First, we note that the component uncertainties does not explain the observed discrepancies. The poor matching between the experimental results and the theoretical prediction is rather imputable to the nonlinearity provided by the nonlinear resistor, which does not follow exactly the sinusoidal law (62). Nevertheless, the multistable property of the system is experimentally established. Indeed, there exist 4 threshold values, U th1 , U th2 , U th3 and U th4 which allow to determine the final state of the elementary cell among the 5 possible stable steady states,

U * 1 , U * 2 , U * 3 , U * 4 , U * 5 .
Therefore, the image processing task inspired by the multistable property of the system is implemented with the electronic device of Fig. 39.

V. CONCLUSION

In this chapter, we have reported a rich variety of image processing operations inspired by the properties of nonlinear systems. Considering a mechanical analogy, we have split the class of nonlinear systems into purely inertial systems and overdamped systems. Using this original description, we have established the properties of nonlinear systems in the context of image processing.

For purely inertial systems, image processing tasks such as contrast enhancement, image inversion, gray level extraction or image encryption can be performed. The applications of the nonlinear techniques presented in this review are similar to those developed by mean of chemical active media [START_REF] Teuscher | Proceedings of the 2005 Workshop on Unconventional Computing, From cellular Automata to Wetware[END_REF], even if these last media are rather overdamped than inertial. In particular, the dynamics of the nonlinear oscillators network, which enables contrast enhancement, can also be used to reveal "hidden images". Indeed, "hidden images" are defined as fragment of picture with brightness very close to the brightness of the image background. Despite a weak difference of brightness between the hidden image and the image background, our nonlinear oscillators network will take advantage of its properties to reveal the hidden image. Another interesting property of this network is to reveal consecutively fields of the image with increasing or decreasing brightness at different processing times. We trust that this feature, also shared by Belousov-Zhabotinsky chemical media, may have potential applications in image analysis in medicine [START_REF] Teuscher | Proceedings of the 2005 Workshop on Unconventional Computing, From cellular Automata to Wetware[END_REF]. Lastly, the noise effects in this purely inertial network leads to cryptography applications. Unlike classical cryptography devices, built with chaotic oscillators, we have proposed an encryption scheme based on the reversibility of our inertial system. Moreover, the encryption key, which ensures the restoration of the initial data, is the time of evolution of the data loaded in the nonlinear network. Therefore, the main advantage of our device is to allow an easy change of the encryption key.

On the other hand, the properties of strongly dissipative or overdamped systems can also give rise to interesting image processing tools. For instance, we have shown the possibility to realize noise filtering, edge detection, or extraction of regions of interest of a weakly contrasted picture. Concerning noise filtering applications based on reaction-diffusion media, the processing is based on the transient behavior of the network since the filtered image depends on the processing times. By contrast, the extraction of regions of interest presents the main advantage to be independent from the processing time since the filtering is realized when the network reaches a stationary pattern. Therefore, this feature can allow an automatic implementation of the processing task.

VI. OUTLOOKS

A. Outlooks on microelectronic implementation

For each nonlinear processing example, we have attempted to propose an electronic implementation using discrete electronic components. Even if these macroscopic realizations are far from real practical applications, they present the main advantage to valid the concept of integration of CNN for future development in microelectronics.

Indeed, in recent years, the market for solid-state image sensors has been experiencing explosive growth due to the increasing demands for mobile imaging systems, video cameras, surveillance, or biometrics. Improvements in this growing digital world continue to be made with two main image sensor technologies: charge coupled devices (CCD) and CMOS sensors. The continuous advances in CMOS technology for processors and DRAMs have made CMOS sensor arrays a viable alternative to the popular CCD sensors. New technologies provide the potential for integrating a significant amount of VLSI electronics into a single chip, greatly reducing the cost, power consumption, and size of the camera [START_REF] Fossum | Active pixel sensors: Are CCDs dinosaurs?[END_REF][START_REF] Fossum | CMOS Image Sensor : Electronic Camera On A CHIP[END_REF][START_REF] Seitz | Solid-State Image Sensing[END_REF][START_REF] Litwiller | CCD vs. CMOS: Facts and Fiction[END_REF]. In past years, most of the work on complex CMOS systems have dealt with the integration of sensors providing a processing unit at chip level (system-on-chip approach) or at column level by integrating an array of processing elements dedicated to one or more columns [START_REF] Yadid-Pecht | In-Pixel Autoexposure CMOS APS[END_REF][START_REF] Acosta-Serafini | A 1/3 VGA Linear Wide Dynamic Range CMOS Image Sensor Implementing a Predictive Multiple Sampling Algorithm With Overlapping Integration Intervals[END_REF][START_REF] Kozlowski | Pixel Noise Suppression via SoC Management of Target Reset in a 1920 × 1080 CMOS Image Sensor[END_REF][START_REF] Sakakibara | A High-Sensitivity CMOS Image Sensor With Gain-Adaptative Column Amplifiers[END_REF]. Indeed, pixel-level processing is generally dismissed because pixel sizes are often too large to be of practical use. However, as CMOS image sensors scale to 0.18 µm processes and under, integrating a processing element at each pixel or group of neighboring pixels becomes feasible. More significantly, employing a processing element per pixel offers the opportunity to achieve massively parallel computations and thus the ability to implement full image systems requiring significant processing such as digital cameras and computational sensors [START_REF] Loinaz | A 200mv 3.3v CMOS Color Camera IC Producing 352 × 288 24-b Video at 30 Frames/s[END_REF][START_REF] Smith | A single-chip 306x244-pixel CMOS NTSC video camera[END_REF][START_REF] Gamal | Pixel level processing -Why, What and How?[END_REF]. The last significant progress in CMOS technologies make possible the realization of vision systems on chip (VSoCs). Such VSoCs are eventually targeted to integrate within a semiconductor substrate the functions of optical sensing, image processing in space and time, highlevel processing, and the control of actuators. These chips consist of arrays of mixed-signal processing elements (PEs) which operate in accordance with single instruction multiple data (SIMD) computing architectures.

The main challenge when designing a SIMD pixel parallel sensor array is the design of a compact, low-power, but versatile and fully programmable processing element. For this purpose, the processing function can be based on the paradigm of Cellular Neural Networks (CNN). CNN can be viewed as a very suitable framework for systematic design of image processing chips [START_REF] Roska | Review of CMOS Implementations of the CNN Universal Machine-Type Visual Microprocessors[END_REF]. The complete programmability of the interconnection strengths, its internal image-memories, and other additional features make this paradigm a powerful front-end for the realization of simple and mediumcomplexity artificial vision tasks [START_REF] Espejo | A 0,8um CMOS programmable analog-array-processing vision-chip with local logic and image-memory[END_REF]. Some proof-of-concept chips operating on preloaded images have been designed [START_REF] Rekeczky | CNN-based spatio-temporal nonlinear ïň Ąltering and endocardial boundary detection in echocardiography[END_REF][START_REF] Czuni | Motion segmentation and tracking with edge relaxation and optimization using fully parallel methods in the cellular nonlinear network architecture[END_REF]. Only a small amount of researches have integrated CNN on real vision chips. As an example, Espejo [START_REF] Espejo | A 64x64 CNN universal chip with analog and digital I/O[END_REF] reports a 64 × 64 pixel programmable computational sensor based on a CNN. This chip is the first fully operational CNN vision-chip reported in literature which combines the capabilities of image-transduction, programmable image-processing and algorithmic control on a common silicon substrate. It has successfully demonstrated operations such as low-pass image filtering, corner and border extraction, and motion detection. More recently, other works have focused on the development of CMOS sensors including the CNN paradigm (Petras, 2003;Carmona, 2003). The chip consists of 1024 processing units arranged into a 32 × 32 grid, and contains approximatively 500000 transistors in a standard 0.5 µm CMOS technology. However, in these pioneering vision chips, the pixel size is often over than 100 µm x 100 µm. Obviously, these dimensions can not be considered as realistic dimensions for a real vision chip. A major part of this crucial problem should be resolved in future years by using the new emergent CMOS technologies. Indeed, CMOS image sensors directly benefit from technology scaling by reducing pixel size, increasing resolution and integrating more analog and digital functionalities on the same chip with the sensor. We expect that further scaling of CMOS image sensor technology and improvement in their imaging performances will eventually allow the implementation of efficient CNN dedicated to nonlinear image processing.

B. Outlooks on processing applications

The nonlinear processing tools developed in this chapter are inherited from the properties of homogeneous media. In the case of applications based on the properties of reaction-diffusion media, it could be interesting to consider the effects of both nonlinearity and structural inhomogeneities. Indeed, novel properties inspired by biological systems, which are rather inhomogeneous than homogeneous [START_REF] Keener | Homogeneization and propagation in the bistable equation[END_REF]Morfu et al, 2002a;Morfu et al, 2002b;[START_REF] Morfu | Propagation failure reduction in a Nagumo chain[END_REF], could allow to optimize the filtering tools developed in this chapter.

For instance, in section IV.C.1, the noise removal method based on the homogeneous Nagumo equation provides a blurry filtered image. It constitutes an additional difficulty to extract the edge of the image with an accurate location. Indeed, noting that the contours of the image corresponds to step-like profiles, the diffusive process increases the spatial expansion of the contours. To avoid this problem, anisotropic diffusion has been introduced in order to reduce the diffusive effect across the image contour. This method has been proposed by Perona and Malik to encourage intraregion smoothing in preference to interregion smoothing [START_REF] Perona | Scale-Space and edge detection using anisotropic diffusion[END_REF]. To obtain this property, Perona and Malik replaced the classical linear isotropic diffusion equation

∂I(x, y, t) ∂t = div(∇I), (68)
by

∂I(x, y, t) ∂t = div(g( ∇I )∇I), (69)
in order to adapt the diffusion with the image gradient. In eqs. ( 68) and (69), I(x, y, t) represents the brightness of the pixel located at the spatial position (x, y) for a processing time t, while ∇I is the gradient amplitude. Moreover, the anisotropy is ensured by the function g( ∇I ) which "stops" the diffusion across the edges. For instance, Perona and Malik considered the function

g(x) = 1 1 + x 2 K 2 , ( 70
)
where K is a positive parameter.

Noting that when x → ∞, g(x) → 0, the effect of anisotropic diffusion is to smooth the original image while the contours are preserved. Indeed, the edge of the image corresponds to brightness discontinuities which lead to strong values of the image gradient [START_REF] Black | Robust anisotropic diffusion[END_REF]. This interesting property of anisotropic diffusion is illustrated in Fig 41 .  For a sake of clarity, the algorithm developed by Perona and Malik is rather extensively detailed in Appendix D and we only discuss here the results obtained by filtering the noisy picture of fig 41 .(a). Contrary to the isotropic nonlinear diffusion based on the Nagumo equation, the edge of the image remains well localized for all the processing times presented in Fig. 41. However, although the noise seems removed for processing times exceeding t = 5, the contrast of the image is never enhanced. Therefore, anisotropic diffusion and nonlinear diffusion do not share the same weakness and it could be interesting to attempt to circumvent the limitations of this two techniques.

For instance, if we compare the continuous equation ( 9) of nonlinear diffusion with the anisotropic equation ( 69) proposed by Perona and Malik, it is clear that the anisotropy can be introduced in our system via the coupling parameter D. Moreover, with the Perona and Malik's method, the pixels brightness does not directly experiences the nonlinearity as in our method. Therefore, the nonlinear noise filtering tool presented in section IV.C.1 could be more efficient if the interesting properties of anisotropic diffusion were also considered by introducing a coupling law. Especially, we expect that the anisotropy preserves the location of the image edges, while the nonlinearity enhances the image contrast and removes the noise in the same time.

Lastly, we close this chapter by presenting another interesting and nonintuitive phenomenon, which takes place in nonlinear systems under certain conditions. This effect, known as Stochastic Resonance effect (SR), has been introduced in the eighties to account for the periodicity of ice ages [START_REF] Benzi | Stochasic resonance in climatic change[END_REF]. Since then, SR has been widely reported in a growing variety of nonlinear systems [START_REF] Gammaitoni | Stochastic Resonance[END_REF], where it has been shown that adding an appropriate amount of noise to a coherent signal at a nonlinear system input enhances the response of the system. Detection of subthreshold signal using noise has been proven in neural information process [START_REF] Longtin | Stochastic resonance in neuron models[END_REF][START_REF] Nozaki | effects of colored noise on stochastic resonance in sensory neurons[END_REF][START_REF] Stocks | Generic noise-enhanced coding in neuronal arrays[END_REF] and in data transmission fields [START_REF] Barbay | Noise-assisted transmission of binary information: theory and experiment[END_REF][START_REF] Zozor | Stochastic resonance in locally optimal detectors[END_REF][START_REF] Morfu | Digital information receiver based on stochastic resonance[END_REF][START_REF] Comte | Stochastic Resonance: another way to retrieve subthreshold digital data[END_REF][START_REF] Duan | Signal detection for frequency-shift keying via shorttime stochastic resonance[END_REF] as well as information transmission in array of such coupled stochastic resonator [START_REF] Lindner | Noise enhanced propagation[END_REF][START_REF] Chapeau-Blondeau | Noise assisted propagation over a nonlinear line of threshold elements[END_REF][START_REF] Báscones | Pulse propagation sustained by noise in arrays of bistable electronic circuits[END_REF][START_REF] Morfu | Propagation failure reduction in a Nagumo chain[END_REF]. Recent studies have also shown that noise can enhance image perception [START_REF] Simonotto | Visual perception of stochastic resonance[END_REF][START_REF] Moss | Stochastic resonance and sensory information processing: a tutorial and review of application[END_REF], autostereograms interpretation [START_REF] Ditzinger | Noise improves threedimensional perception: Stochastic resonance and other impacts of noise to the perception of autostereograms[END_REF], visual perception via eyes micro-saccade in retina (Hongler et al, 1998), image processing [START_REF] Vaudelle | Stochastic resonance and noise-enhanced transmission of spatial signals in optics: The case of scattering[END_REF][START_REF] Chapeau-Blondeau | Stochastic resonance and the benefit of noise in nonlinear systems, Noise, Oscillators and Algebraic Randomness-From Noise in Communication Systems to Number Theory[END_REF][START_REF] Histace | Constructive action of noise for impulsive noise removal in scalar images[END_REF][START_REF] Blanchard | Constructive action of the speckle noise in a coherent imaging system[END_REF]. The investigation of noise effects in nonlinear systems is undoubtedly of great interest in nonlinear signal processing or in image processing context [START_REF] Zozor | Stochastic resonance in discrete time nonlinear AR(1) models[END_REF][START_REF] Zozor | Noise aidded processing: revisiting dithering in a Sigma-Delta quantizer[END_REF]. We thus propose to present the phenomenon of S.R. using the methodology exposed in [START_REF] Chapeau-Blondeau | Stochastic resonance and the benefit of noise in nonlinear systems, Noise, Oscillators and Algebraic Randomness-From Noise in Communication Systems to Number Theory[END_REF]. Moreover, to show a visual perception of the S.R. effect, we consider the black and white image of fig.42.(a), where we note the probability p 1 to have a white pixel and p 0 = 1 -p 1 the probability to have a black one. A gaussian white spatial noise η i,j with R.M.S. amplitude value σ is added in each pixel I i,j of the initial image. The resulting noisy image is then threshold filtered with a threshold V th to obtain the image I b , according to the following threshold filtering rule: if I i,j + η i,j > V th then I bi,j = 1 else I bi,j = 0.

(71)

The similarity between the two images I and I b can then be quantified by the cross-covariance (Chapeau-Blondeau, 2000)

C or by

R II b = II b I 2 I 2 b , (73)
where < . > corresponds to an average over the images. These two similarity measures are defined by

R II b = p 1 (1 -F η (V th -1)) p 1 p 1 1 -F η (V th -1) + 1 -p 1 1 -F η (V th )
.

and

C II b = p 1 (1 -F η (V th -1)) -p 1 q 1 (p 1 -p 2 1 )(q 1 -q 2 1 )
.

with q 1 = p 1 (1 -F η (V th -1)) + (1 -p 1 )(1 -F η (V th )) and where F η is the cumulative distribution function of the noise.

In the case of a gaussian white noise of R.M.S. amplitude σ, the cumulative distribution function can be expressed under the form

F η (u) = 1 2 + 1 2 erf u √ 2σ . (74)
In eq. ( 74) the error function is defined by erf (u) = 2 √ π u 0 exp(-t 2 )dt. The two quantities (72) and ( 73 Even if the model of human visual perception is more complex than a standard threshold filtering [START_REF] Bálya | CN N framework for modeling parallel processing in a mammalian retina[END_REF], this simple representation is convenient to determine analytically the optimum amount of noise that provides the best visual perception of images via Stochastic Resonance. Moreover, the S.R. phenomenon is shared by a wide class of nonlinear systems including neural networks which also intervenes in the process of images perception. Since neurons are basically threshold devices that are supposed to work in a noisy environment, the interest of considering noise effect seems to be of crucial importance in order to develop artificial intelligence applications that perfectly mimic the real behavior of nature. Therefore, for the next decades, we trust that one of the most interesting challenge could be to complete the description of nonlinear models by including the contribution of noise effects.

Assuming that when t → +∞, the particle evolves to the steady states W = 0 for W 0 < 1/2 and W = 1 for W 0 > 1/2, we finally obtain the displacement of the particle with initial position W 0 as

W (t) = 1 2 1 + W 0 -1 2 (W 0 -1 2 ) 2 -W 0 (W 0 -1)e -1 2 t . ( 82
)

APPENDIX B. RECALL ON JACOBIAN ELLIPTIC FUNCTION

We recall here the properties of Jacobian elliptic functions used in section II.B. These three basic functions cn(u, k), sn(u, k), dn(u, k) play an important role in nonlinear evolution equations and arise from the inversion of the elliptic integral of first kind [START_REF] Abramowitz | Handbook of Mathematical Functions[END_REF]:

u(ψ, k) = ψ 0 dz 1 -k sin 2 z , (83)
where k ∈

Considering the circular function properties, we also have

sn 2 (u, k) + cn 2 (u, k) = 1. ( 86
)
Moreover, using the result (84), we obtain the following identity:

dn 2 (u, k) + ksn 2 (u, k) = 1. ( 87
)
where W 0 denotes the initial position of the particle.

Inverting the tangent function, provides straightforwardly the solution of eq. ( 88) under the form

W (t) = 1 π(n -1)
arctan tan(π(n -1)W 0 )e -β(n-1)

2 2πt + k n -1 , ( 94
)
where k is an integer coming from the tangent inversion.

Note that from a physical point of view, k must ensure that the particle position evolves towards one of the stable states of the system for a sufficiently large time, that is when t → +∞. Indeed, for an initial condition between two consecutive unstable steady states, the asymptotic behavior of the uncoupled network can reduce to the following rule:

if 2k -1 2(n -1) < W 0 < 2k + 1 2(n -1) W (t → +∞) = k (n -1) (95)
This rule can be transformed to get

if k - 1 2 < (n -1)W 0 < k + 1 2 W (t → +∞) = k (n -1) (96)
Lastly, identifying eq. ( 94) with eq. ( 96) when t → +∞, we deduce that k must be the nearest integer of W 0 (n -1).

APPENDIX D. PERONA AND MALIK ANISOTROPIC DIFFUSION

ALGORITHM

We recall here the algorithm introduced by Perona and Malik to compute their method based on anisotropic diffusion equation. The anisotropic diffusion equation ( 69) can be discretized with the time step dt to obtain

In eq. ( 97), I t s represents the brightness of the pixel located at the position s in a discrete two-dimensional grid which corresponds to the filtered image after a processing time t. η s is the number of neighbors of the pixel s, that is 4, except for the image edge where η s = 3 and for the image corners where η s = 2. The spatial neighborhood of the pixel s is noted N r. The local gradient ∇I s,p can be estimated by the difference of brightness between the considered pixel s and its neighbor p: ∇I s,p = I p -I t s , p ∈ N r.

(98)

Lastly, the description of the system is completed by defining the edge stopping function g(x) as the lorentzian function:

g(x) = 1 1 + x 2 K 2 , ( 99
)
where K is a positive parameter.

FIGURE 2 .

 2 FIGURE 2. Bistable behaviour of the overdamped system in the case α = 1/2. Left: Evolution of a particle for different initial conditions in the range [0; 1]. The solid line is plotted with the analytical expression (8) whereas the (o) signs correspond to the numerical solution of eq. (5) for different initial conditions W 0 ∈ [0; 1]. The potential φ obtained by integrating the nonlinear force (3) is represented at the right to provide a reference.

  a kink propagates from left to right with a positive velocity c = D/2(1 -2α) (fig. 4.(a) left). 2. Else, if the profile is set by W (ξ-ξ 0 ) , a kink propagates from right to left with a negative velocity c = -D/2(1-2α) (fig. 4.(a) right).

  a kink propagates from right to left with a negative velocity c = D/2(1 -2α) (fig. 4.(b) left).

FIGURE 4 .

 4 FIGURE 4. Propagative solution of the continuous Nagumo equation (9) with D = 1. Spatial representation of the kink for t = 0 in dotted line and for t = 20 in solid line. The arrow indicates the propagation direction, the corresponding potential is represented at the right end to provide a reference. (a) α = 0.3. (b) α = 0.7.

  a kink propagates from left to right with a positive velocity c = -D/2(1 -2α) (fig. 4.(b) right).

FIGURE 5 .

 5 FIGURE 5. Double well potential deduced from the nonlinear force (11) represented for m = 2.58, α = 1.02 and ω 0 = 1. A particle with an initial condition W 0 i < mα √ 2 evolves with an initial potential energy above the barrier ∆.

as shown in figure 6 FIGURE 6 .

 66 FIGURE 6. . (a): Normalized pulsation ω/ω 0 versus the initial condition W 0 i . (b) Modulus parameter k versus W 0 i . The parameters of the nonlinearity m = 2.58, α = 1.02 impose the allowed amplitude range ] -∞; 1.137] [4.023; + ∞[.

FIGURE 7 .

 7 FIGURE 7. (a) Temporal evolution of the two oscillators. Top: evolution of the first oscillator with initial condition W 0 1 = 0. Bottom: evolution of the second oscillator with initial condition W 0 2 = 0.2. (b) Temporal evolution of the displacement difference δ between the two oscillators. Parameters: m = 2.58, α = 1.02 and ω 0 = 1.

FIGURE 8 .

 8 FIGURE 8. Influence of the nonlinearity parameter α on the displacement difference δ between the two oscillators of respective initial conditions 0 and 0.2. Parameters m = 2.58 and ω 0 = 1. (a): (topt = 1.75 × 10 -3 ; α = 0.4). (b): (topt = 1.66 × 10 -3 ; α = 1.05). (c): (topt = 1.25 × 10 -3 ; α = 1.5). (d): (topt = 0.95 × 10 -3 ; α = 1.63).

FIGURE 9 .

 9 FIGURE 9. Filtered images and their corresponding histogram obtained with the nonlinear oscillators network (23) for different processing times. (a) : Initial image (t = 0). (b) : t = 3.96 × 10 -5 . (c) : t = 1.4 × 10 -4 . (d) : t = 1.81 × 10 -4 . (e) : t = 3 × 10 -4 . (f ) : t = 3.5 × 10 -4 . (g) : t = 4.64 × 10 -4 . (h) : t = 5.21 × 10 -4 . (i) : t = 6.29 × 10 -4 . (j) : t = 1.47 × 10 -3 . (k) : t = topt = 1.64 × 10 -3 . Parameters: m = 2.58, α = 1.02, ω 0 = 1.

  FIGURE 10. Response of the nonlinear system for different nonlinearity parameters α at the corresponding optimal time topt (solide line) compared to a uniform rescaling (dotted line). The curves are obtained with eqs. (20) and (21) setting the time to the optimum value defined by the maximum of δ(t) (see fig.8). In addition, we let the initial conditions W 0 i varying in the range [0; 0.2] in eqs. (20) and (21). (a): (topt = 1.75 × 10 -3 ; α = 0.4). (b): (topt = 1.66 × 10 -3 ; α = 1.05). (c): (topt = 1.25 × 10 -3 ; α = 1.5). (d):(topt = 0.95 × 10 -3 ; α = 1.63). ω 0 = 1

FIGURE 11 .

 11 FIGURE 11. Gray level extraction. The response of the system is represented at the top of each figure. At the bottom of each figure, a threshold filtering of the filtered image is realized replacing the pixel gray level with 0.2 (white) if that gray level exceeds the threshold V th = 0.25, otherwise with 0 (black). (a) : Initial gray scale (t = 0). (b) : t = 3.3 × 10 -3 . (c) : t = 3.33 × 10 -3 . (d) : t = 3.36 × 10 -3 . (e) : t = 3.39 × 10 -3 . (f ) : t = 3.42 × 10 -3 . (g) : t = 3.45 × 10 -3 . (h) : t = 3.48 × 10 -3 . (i) : t = 3.51 × 10 -3 . (j) : t = 3.54 × 10 -3 . Nonlinearity parameters: m = 2.58, α = 1.02 and ω 0 = 1

  FIGURE 12. Encrypted image and the corresponding response of the nonlinear oscillators network for different times exceeding topt. (a) : Enhancement of contrast of the initial image for t = topt = 1.64 × 10 -3 . (b) : t = 3.28 × 10 -3 . (c) : t = 6.56 × 10 -3 . (d) : t = 8.24 × 10 -3 . (e) : t = 9.84 × 10 -3 . Parameters: m = 2.58, α = 1.02, ω 0 = 1.

FIGURE 14 .

 14 FIGURE 14. Theoretical cubic law (24) in solid line compared to the experimental characteristic plotted with crosses. Parameters: m = 2.58 V , α = 1.02 V , K = 10.

  FIGURE 15. : (a): Temporal evolution of two elementary cells of the chain with respective initial conditions W 0 1 = 0V (top) and W 0 2 = 0.2V (bottom). (b): Evolution of the voltage difference between the two oscillators. Parameters: K = 10, R = 10KΩ, C = 10nF , m = 2.58V , α = 1.02V , topt = 1.46ms.

  . (a) for two different initial conditions, namely W 0 1 = 0V (top) and W 0 2 = 0.2V (bottom). As previously explained, the way to introduce the initial condition allows to balance the dissipative

  FIGURE 16. Response of the system to a set of initial conditions W 0 i ∈ [0; 0.2] at the optimal time. The solid line is obtained with eqs. (20), (21) and (27) setting the time to the theoretical optimal value 1.64ms, the initial condition varying in [0; 0.2V ]. The crosses are obtained experimentally for the corresponding optimal time 1.46ms. Parameters: R = 10KΩ, C = 10nF , m = 2.58V , α = 1.02V , K = 10.

FIGURE 17 .

 17 FIGURE 17. Theoretical response of the purely inertial system (solid line) compared to the experimental ones (crosses) for 4 different times and for a range of initial conditions [0; 0.5V ]. Parameters: R = 10KΩ, C = 10nF , m = 2.58V , α = 1.02V , K = 10. (a) experimental time t = 564µs corresponding to the theoretical time t = 610µs. (b) experimental time t = 610µs and theoretical time 713µs. (c) experimental time t = 675µs and theoretical time 789µs. (d) experimental time t = 720µs and theoretical time 841µs.

  FIGURE 18. (a) Nonlinear electrical lattice. (b): The nonlinear resistor R NL .

FIGURE 20 .

 20 FIGURE 20. Realization of a nonlinear resistor with a polynomial generation circuit. β = 10VaV b

  FIGURE 21. . Current-Voltage characteristics of the nonlinear resistor of fig. 20. Parameters: β = -10VaV b , Va = -2V , V b = 2V .

  FIGURE 22. (a) Temporal evolution of a uniform initial condition U 0 = 0.4 applied to the whole network. (b) Temporal evolution of the perturbation applied to the cell n = 24 for D = 0.5 and b 0 = 0.2. (c) Temporal evolution of the perturbation applied to the cell n = 24 for D = 5 and b 0 = 0.2. Solid line: theoretical expressions (41) and (51), (o) signs: numerical results.

  FIGURE 23. Response of the lattice to a uniform initial condition corrupted by a constant perturbation at two different processing times. (o) signs: numerical results; solid line: theoretical expression (51). (a): D = 5; (b): D = 0.5. (I) initial condition for t = 0, (II) state of the lattice for t = 1 , (III) state of the lattice for t = 2.

  FIGURE 24. Spatio-temporal view of the response of the lattice to the previous initial condition. (a): D = 5. (b): D = 0.5.

FIGURE 26 .

 26 FIGURE 26. (a) Temporal evolution in normalized units of a uniform initial condition W 0 = 0.4 applied to the network. (b) Temporal evolution of the perturbation applied to the cell n = 24 for b 0 = 0.2 and D = 0.5 corresponding to a coupling resistor R = 3KΩ. (c) Temporal evolution of the perturbation applied to the cell n = 24 for b 0 = 0.2 and D = 5 corresponding to a coupling resistor R = 300Ω. C = 33nF , Nonlinearity parameters β = 1, V b = 1.12V , Va = 0.545V involving α = 0.49.

  FIGURE 27. Response of the lattice to a uniform initial condition corrupted by a constant perturbation at two different processing times. Parameters: C = 33nF , V b = 1.12V , Va = 0.545V , α = 0.49. (a): R = 300Ω that is D = 5;(b): R = 3KΩ that is D = 0.5. (I) initial condition for t=0, (II) state of the lattice for t = 2 (τ = 0.1ms), (III) state of the lattice for t = 4 (τ = 0.2ms).

  FIGURE 28. Noise filtering of a one dimensional signal with an electrical nonlinear lattice. (a) : normalized noisy sinusoidal signal given by eq. (52) loaded as initial condition at the nodes of the lattice. σ = 0.15, N = 48 and A = 0.264. (b) : Filtered signal obtained for a processing time t = 0.6 (τ = 92.3µs) and a coupling D = 0.5 (that is R = 3KΩ). (c) : Filtered signal obtained for a processing time t = 0.6 (τ = 92.3µs) and a coupling D = 5 (that is R = 300Ω). Parameters: C = 100nF , β = 1, V b = 1.12V , Va = 0.545V .

FIGURE 29 .

 29 FIGURE 29. Noisy image of the Coliseum.

FIGURE 30 .

 30 FIGURE 30. : Noise filtering of the image represented in Fig. 29. (a), (b), (c), (d). Filtered image obtained for D = 0.075 and for the respective processing times t = 1, t = 3, t = 6 and t = 9. (e), (f ), (g), (h) Filtered image obtained for D = 0.1 and for the respective processing times t = 1, t = 3, t = 6 and t = 9. (i), (j), (k), (l) Filtered image obtained for D = 0.2 and for the respective processing times t = 1, t = 3, t = 6 and t = 9. (m), (n), (o), (p) Filtered image obtained for D = 0.3 and for the respective processing times t = 1, t = 3, t = 6 and t = 9.

  FIGURE 31. Contour detection of a black square in a white background. (a) initial image and its profile. (b) Edge detection of the object and its profile obtained with the standard cubic nonlinearity (5) with threshold α = 1/3. Processing time t = 4, D = 1. (c) contour and the corresponding profile obtained with the nonlinearity (54). Processing time t = 4, D = 1.

  FIGURE 33. Noisy and weakly contrasted image of a soldering between two rods of metal. The image histogram is represented at the right.

W 0 fFIGURE 34 .

 034 FIGURE 34. Mechanical point of view of the bistable overdamped network used for image processing. (a): The pixel with coordinates i, j and gray level W i,j is analog to an overdamped particle coupled by springs of strength D to its 4 nearest neighbors. (b) : The particle is attracted in one of the two wells of the bistable potential according to the resulting elastic force applied by the 4 coupled particles.

FIGURE 36 .

 36 FIGURE 36. Multistable potential represented for β = 9.82 × 10 -2 and n = 5. The potential barrier between two consecutive extrema is β/π.

FIGURE 37 .

 37 FIGURE 37. Temporal evolution of an overdamped particle experiencing the multistable potential. Parameters: n = 5 and β = 0.25. Solid line: theoretical expression (60), (o) signs: numerical results obtained solving eq. (57)

  Fig. 40.(b) to the theoretical expression (62). The weak discrepancies observed between the theoretical and experimental laws can be reduced by increasing the order

FIGURE 40 .

 40 FIGURE 40. (a) : Response of an elementary cell of the multistable network to different initial conditions in the uncoupled case. The plot of the theoretical expression (60) in solid line is compared to the experimental results represented by crosses. (b) : Nonlinear current-voltage characteristics. The sinusoidal law (62) in solid line matches the experimental characteristics in + signs. The Components values are R 0 = 2KΩ, C = 390nF , I M = 2mA. The zeros of the sinusoidal current defines the 4 unstable states U th1 , U th2 , U th3 , U th4 which correspond to thresholds and the 5 stable steady states U * 1 , U * 2 , U * 3 , U * 4 , U * 5 which correspond to attractors.

  FIGURE 41. : Noise filtering based on anisotropic diffusion. The filtering image are obtained using the algorithm detailed in Appendix D with the parameters dt = 0.01 and K = 0.09. (a) initial image, (b), (c), (d), (e), (f ), (g), (h), (i) image for the respective processing times t = 1, t = 2, t = 3, t = 4, t = 5, t = 6, t = 7, t = 8.

  FIGURE 42. : (a) : Initial black and white image with p1 = 0.437. (b) : similarity measures (72) and (73) versus the noise RMS amplitude value σ for V th = 1.1.

  ) are plotted versus the RMS noise amplitude σ in figure 42.(b), where a resonant-like behavior reveals the standard stochastic resonance signature. Indeed, there exists an optimum amount of noise that maximizes the similarity measures (72) and (73). According to fig.42.(b), this optimal noise RM S value is σ = 0.4. To valid the similarity measures, we qualitatively analyze the pictures obtained for different noise amplitudes . It is confirmed in fig.43that the noise optimal value σ = 0.4 allows the best visual perception of the coliseum through the nonlinear systems.

  FIGURE 43. :(a), (b), (c), (d), Threshold filtered image with the rules (71) and a threshold V th = 1.1 and with white gaussian noise with respective RMS noise amplitude σ = 0.1, σ = 0.4, σ = 0.8, σ = 1.4.

  [0; 1] is the elliptic modulus. The jacobian elliptic functions are defined by sn(u, k) = sin(ψ), cn(u, k) = cos(ψ), dn(u, k) = 1 -k sin 2 (ψ). (84) This definition involves the following properties for the derivatives: d sn(u, k) du = cn(u, k)dn(u, k), d cn(u, k) du = -sn(u, k)dn(u, k), d dn(u, k) du = -ksn(u, k)cn(u, k),

  s,p )∇I s,p .

ACKNOWLEDGEMENTS

First of all, S. Morfu would like to thank J.M. Bilbault and O. Michel for giving him the opportunity to evolve in this fascinating scientific world. S. Morfu would like to extend his appreciation to P.O. Amblard, F. Chapeau-Blondeau and D. Rousseau who have brought specific references to his attention. S. Morfu is also grateful to his colleague Julien Dubois for useful advices.

P. Marquié would like to dedicate this chapter to Arnaud and Julie and to thanks J.M. Bilbault and M. Remoissenet.

The authors would like to take the opportunity to highlight the technical assistance of M. Rossé for the design of the sinusoidal nonlinear resistor.

Lastly, S. Morfu would like to warmly dedicate this chapter to Giovanni and Grazia Morfu.

APPENDIX A. RESPONSE OF A CELL OF THE OVERDAMPED NETWORK

In the uncoupled case, and for α = 1/2, a particle of displacement W obeys to

Separating the variables in Eq. ( 75) gives

which can be integrated to obtain

where K is an integration constant. Equation ( 77) can be arranged to become a second order equation in W

Provided that the discriminant is positive, the solutions are given by:

Assuming that initially, the position of the particle is W (t = 0) = W 0 , we can express the integration constant K under the form

Inserting the constant (80) in the solution (79) leads to the following expression of the displacement:

The equation of motion of an overdamped particle submitted to the sinusoidal force (56) expresses

where W represents the particle displacement.

The steady states of the system are deduced from the zeros of the nonlinear force. Using the methodology exposed in section II.A.1, we can establish that the roots of the nonlinear force correspond alternatively to unstable and stable steady states. If k is an integer, the unstable and stable states of the system respectively write

Separating the variables of eq. ( 88), we obtain

Using the identity sin(2a) = 2 sin a cos a, eq. ( 90) becomes

Next, considering the derivative of the tangent function in eq. ( 91), yields

A direct integration of eq. ( 92) gives tan π(n -1)W = tan π(n -1)W 0 e -β(n-1) 2 2πt .

(93)