Arthur Robinet
email: a.robinet@brgm.fr

Bruno Castelle

Déborah Idier

Maurizio D'anna

Gonéri Le Cozannet

Simulating the Impact of Sea-level Rise and Offshore Bathymetry on Embayment Shoreline Changes

Keywords: Shoreline model, sediment transport, bathymetry, waves, sea-level rise

LX-Shore is a reduced-complexity shoreline change model driven by cross-shore and longshore processes which can account for man-made or natural non-erodible areas such as groynes and headlands. Here we describe and further test the implementation of two recent developments allowing to account for (i) real and non-erodible offshore bathymetric features such as rocky outcrops and canyons affecting offshore wave transformation and, in turn, shoreline variability and (ii) shoreline retreat due sea-level rise. After a description of the numerical developments, the benefits of these new developments are demonstrated with the application of LX-Shore to an idealized embayed beach exposed to real wave climate during a 10-yr period. Three simulations are conducted to test the impact of an outcrop in the middle of the embayment and of a gradual 1-m sea-level rise on shoreline spatial and temporal modes of variability. Results show that the equilibrium planview shoreline and the shoreline variability are strongly impacted by only slightly modifying the bathymetry and varying the mean sea level. These results show the potential of LX-Shore to better understand and further predict shoreline change along real coasts exhibiting mixed (sandy/rocky) and complex seabed morphologies and undergoing sea-level rise.

INTRODUCTION

Worldwide, wave-dominated sandy coasts are subject to complex changes driven by a myriad of physical processes acting on a wide range of spatiotemporal scales [START_REF] Vitousek | A model integrating longshore and cross-shore processes for predicting long-term shoreline response to climate change[END_REF]. Although most of the physical processes are well characterized, accurately simulating shoreline changes along these coasts is still challenging. In the last decades, the development of a new generation of shoreline change models allowing predictions of past and future shoreline changes has received a growing interest from the scientific coastal community and coastal risk managers [START_REF] Montaño | Blind testing of shoreline evolution models[END_REF]. Physical-based reduced-complexity shoreline models focus on shoreline changes only and rely on the few processes that are dominant on the timescales from storm event to decades up to century. They are computationally cheap and robust.

However, reduced-complexity models usually oversimplify the bathymetry and/or the wave transformation processes. They cannot therefore account for the effect of complex seabed morphologies, sometimes observed along real coasts, that may impact the characteristics of onshore-propagating waves and, in turn, the nearshore sediment transports and shoreline changes. Recently [START_REF] Limber | Modeling large-scale shoreline change caused by complex bathymetry in low-angle wave climates[END_REF] coupled the Coastal Evolution Model (Ashton and Murray., 2006) with the spectral wave model SWAN (Booij et al., 1999) and included nearshore shoals in the bathymetry. They thus highlighted the control of the seabed morphology on the equilibrium planview shoreline. However, their model only takes into account shoreline change driven by gradients in longshore transport.

Accounting for the sea-level rise (SLR) in reducedcomplexity models is also crucial to tackle future shoreline changes on the timescales of decades and more (e.g. [START_REF] Vitousek | A model integrating longshore and cross-shore processes for predicting long-term shoreline response to climate change[END_REF]. Indeed, it is expected that sea-level rise will dominate shoreline erosion by the end of the century (Le [START_REF] Cozannet | Uncertainties in Sandy Shorelines Evolution under the Bruun Rule Assumption[END_REF]. SLR may also affect the shoreline variability on shorter timescales in presence of complex offshore bathymetries, with increased depths potentially affecting offshore wave propagation, and, in turn, breaking wave conditions and sediment transport along the coast. This contribution presents the numerical developments performed to integrate the above-mentioned processes in the state-of-the-art shoreline change model LX-Shore [START_REF] Robinet | A reduced-complexity shoreline change model combining longshore and cross-shore processes: The LX-Shore model[END_REF]. The benefits of these new developments are illustrated by applying LX-Shore to an idealized embayed beach exposed to a real wave climate during a 10-yr period.

METHODS

Overview of LX-Shore

LX-Shore is a one-line-based reduced-complexity model resolving shoreline change driven by gradients in longshore sediment transport and by cross-shore processes resulting from disequilibrium between the instantaneous wave conditions and the recent past wave conditions to which the beach was exposed. The total longshore sediment transport is computed with the formula of [START_REF] Kamphuis | Alongshore sediment transport rate[END_REF] and the cross-shore processes are described by an adaptation of the ShoreFor model (Splinter www.cerf-jcr.org _________________________________________________________________________________________________ Journal of Coastal Research, Special Issue No. 95, 2020 et al., 2014). LX-Shore does not resolve explicitly the displacements of specific shoreline points, but instead relies on a 2D grid of cells containing sediment fraction values (ratio of dry to wet area in each grid cell) that are updated at the end of every simulation timestep (based on above-mentioned processes). The shoreline is then retrieved using an interface reconstruction method applied to this grid. LX-Shore is dynamically coupled with SWAN to ensure realistic and robust prediction of breaking wave conditions, which are then used to drive longshore and cross-shore processes. In the former version of LX-Shore [START_REF] Robinet | A reduced-complexity shoreline change model combining longshore and cross-shore processes: The LX-Shore model[END_REF], sea-level rise was not taken into account and the bathymetry provided to SWAN was idealized and computed by LX-Shore at every simulation timestep using the approach of Kaergaard and Fresoe (2013), i.e. by imposing a stationary Dean-like beach profile seaward the simulated planview shoreline. Below we describe two recent developments.

Inclusion of Real Seabed Morphologies

One new development is the possibility to include real bathymetric features and to specify non-erodible bottoms in the bathymetry used by SWAN. If a measured bathymetry, hereafter referred to as the reference bathymetry (h r ), is provided to LX-Shore, it will be automatically merged with the idealized bathymetry (h i ) originally computed by LX-Shore at every simulation timestep 1). This processing is based on two threshold depths, the depth of closure (D c ) and a tunable depth (D s ) greater than D c , and on the following rules: (i) h i is kept near the shoreline, i.e. for 0 ≤ h i < D c (Figure 1

h m = h i + (h r -h i ) c ss (1) c ss = 3 NDR 2 -2 NDR 3 (2) NDR = (h i -D c ) / (D s -D c ) (3)
where c ss a sigmoid-shaped coefficient (function of h i , Figure 1b) and NDR a normalized depth ratio ranging from 0 to 1. Lastly, the resulting bathymetry is compared to the reference bathymetry where the bottom is non-erodible (Figure 1, purple line). Adjustments can then be made to ensure that the final bathymetry never gets deeper than the bathymetry of nonerodible bottoms.

Inclusion of the Sea-Level Rise

A sea-level-change module has been implemented in LX-Shore to allow: (i) variations of the mean sea-level above the reference bathymetry initially provided to LX-Shore; (ii) shoreline retreat driven by SLR as predicted by the so-called Bruun Rule (Bruun, 1962). The first approach allows exploring how the shoreline dynamic is affected by changes in the characteristics of incoming waves due to variations of the mean sea-level. With the second approach, the planview shoreline will progressively retreat as the mean sea-level rises based on the following equation:

dS = dη / tan α (4)
where dS is the cross-shore shoreline retreat associated with the sea-level rise dη and tan α is the average beach profile slope calculated between the shoreline position and D c . Two important assumptions are made: (i) the beach profile used to compute the idealized bathymetry remains the same; (ii) the mean-sea-level contour of emerged rocky structures remains the same as these structures are considered vertical and infinitely high. Test Case LX-Shore is applied to an idealized embayed beach (Figure 2a) whose characteristics and wave climate are inspired by those of Narrabeen beach (Sydney, Australia) [START_REF] Turner | A multidecade dataset of monthly beach profile surveys and inshore wave forcing at Narrabeen[END_REF]. The model is forced with the wave conditions (H s , T p , D p ) hindcasted at the Sydney buoy (Durrant et al., 2014) between 22/07/2005 and 22/07/2015. The tide is neglected and the shoreline position is assumed to correspond to the intersection between the mean sea level and the topography. The beach profile used to compute the idealized bathymetry was calibrated using the time and spatial average of bathymetric profiles measured at Narrabeen along cross-shore transects (transects PF2, PF4 and PF6 in [START_REF] Turner | A multidecade dataset of monthly beach profile surveys and inshore wave forcing at Narrabeen[END_REF] down to nearly 14-m depth. Depth of closure was roughly estimated to 8 m using the formula of Hallermeier (1981) based on SWAN-simulated breaking wave conditions extracted at the center of the embayment. No hypothesis was made on the planview shoreline used to start the simulations. To obtain an initial shoreline in relative equilibrium with the wave climate a preliminary simulation (s0) is performed with LX-Shore using a straight N-S-oriented shoreline (Figure 2a, dotted black line) and enabling only the longshore sediment transport. The shoreline rapidly curves and rotates around a pseudo-equilibrium planview shoreline. The latter is taken as the initial condition of each simulations presented in the next section and is computed as the temporal average of the planview shorelines simulated from 2010-2015 in s0 (Figure 2a, red line).

_________________________________________________________________________________________________

Journal of Coastal Research, Special Issue No. 95, 2020

Effects of Submerged Morphologies and Sea-Level Rise

Three simulations are conducted with different model configurations (Table 1) to illustrate the benefits of the new developments made in LX-Shore. In all simulations, both longshore and cross-shore processes are enabled. Breaking wave conditions extracted in the center of the embayment during simulation s0 and the cross-shore shoreline positions measured at the Narrabeen PF6 transect [START_REF] Turner | A multidecade dataset of monthly beach profile surveys and inshore wave forcing at Narrabeen[END_REF] from 2005 to 2015 are used to calibrate the free parameters of the cross-shore model using the methodology presented in [START_REF] Robinet | A reduced-complexity shoreline change model combining longshore and cross-shore processes: The LX-Shore model[END_REF]. In simulation s1, the reference bathymetry provided to LX-Shore is the idealized bathymetry that LX-Shore computes at the first timestep of the simulation (Figure 3a). Simulation s2 is similar to s1 but a Gaussian-shaped outcrop is added to the reference bathymetry in the center of the embayment (Fig 3b). This outcrop culminates in 6-m depth, has a diameter of nearly 500 m and lies in depths ranging from 12 to 25 m. Simulation s3 is similar to s2 but a constant SLR rate is applied throughout the simulation, which leads to a +1-m elevation of the mean sea level at the end of the simulation.

RESULTS

Average Planview Shoreline

Figure 4 shows the initial planview shoreline (dotted grey line) and the averages of the planview shorelines simulated in s1 (thick black line), s2 (blue line) and s3 (red line). The overall shape of the averaged planview shorelines is similar, though substantial changes are visible between y = 0 m and y = 2000 m. The inclusion of the outcrop in simulation s2 led to an adjustment of the averaged planview shoreline with a mean shoreline retreat (progradation) of about 15 m for 1300 < y < 2000 m (y ≤ 1300 m) and a maximum retreat of 30 m near y = 1700 m. Although not clearly visible, the inclusion of the SLR effects in simulation s3 led to a shoreline retreat of about 10 m all along the beach.

Modes of Shoreline Variability

Following the methodology used in Harley et al. (2015), an Empirical Orthogonal Function (EOF) decomposition is applied to the demeaned timestacks of simulated planview shorelines to extract the three first spatiotemporal modes of shoreline variability simulated in s1, s2, and s3 (thick black lines, blue lines and red lines in Figure 5, respectively).

Interestingly, for simulation s1 the two first dominant modes are similar to those identified at Narrabeen in Harley et al.

_________________________________________________________________________________________________

Journal of Coastal Research, Special Issue No. 95, 2020 (2015). They are the beach rotation mode (black line in Figure 5a,d) related to longshore transport, and the beach translation mode (black line in Figure 5b,e) related to cross-shore transport and whose spatial magnitude decrease toward the south. The third mode, however, is different and seems here to correspond to a process of beach curvature variation resulting from transfer of sediments contained between 1500 < y < 2700 m to adjacent extremities on the timescales from storm event to years. The inclusion of the outcrop in simulation s2 induces significant changes mostly in the spatial component of the modes of shoreline variability (Figure 5a,b,c) essentially located near the center of the beach. One major change with respect to simulation s1 is the apparent attenuation of the translation process near yy = 1800 m (Figure 5b). The presence of this outcrop also appears to slightly inhibit (increase) the contribution of the rotation (translation) mode to overall shoreline variability (Figure 5a,b).

The inclusion of the SLR effects in simulation s3 modifies essentially the spatiotemporal components of mode 2 and 3 (Figure 5b,c,e,f). The attenuation of the translation mode simulated in simulation s2 is removed (Figure 5b) and the trend in the temporal component of mode 2 is reversed with respect to simulations s1 and s2 (Figure 5e).

DISCUSSION

Consistently with the philosophy of reduced-complexity models, simple and computationally cheap developments have been implemented in LX-Shore to account for the complex bathymetry and sea-level rise effect. Here, we presented a test case derived from the Narrabeen beach. A great similarity is found between the dominant modes of shoreline variability simulated by LX-Shore (s1 case) and those isolated at Narrabeen (Harley et al., 2015). Thus, despite its apparent simplicity, LX-Shore appears as a relevant tool to simulate shoreline dynamics and to explore the respective contribution of the physical processes involved in shoreline changes.

The comparative analysis of the results of simulations s1 and s2 suggests that the presence of the Gaussian-shaped outcrop at the center of the embayment modifies the average planview shoreline, which is consistent with the findings of [START_REF] Limber | Modeling large-scale shoreline change caused by complex bathymetry in low-angle wave climates[END_REF]. The outcrop also appears to strongly impact the shoreline spatial and temporal variability, as evidenced by the changes in the spatial component of the dominant modes of shoreline variability. Although the shoreline dynamic is predominantly modified in the lee of the outcrop from a modalwave-direction point of view (Figure 4 and 5), non-negligible adjustments are also induced elsewhere along the embayment (e.g. Figure 5b,c). These results suggest that the sea bottom morphology, and not only the combination of the wave climate and the headland sheltering, contributes to the shoreline dynamic along embayed beaches presenting complex bathymetries (e.g. presence of outcrops and/or submerged rocky platforms typically extending from headlands).

The sea-level rise also has a significant effect on the simulated shoreline variability. While the beach rotation mode appears to be insensitive to this process, this is not the case of the translation mode whose spatial component is similar to the spatial component of the translation mode obtained when the outcrop is absent from the bathymetry. The continuous increase in the mean

_________________________________________________________________________________________________

Journal of Coastal Research, Special Issue No. 95, 2020 sea level during the simulation may have damped the wave transformation processes occurring above the outcrop, resulting in more uniform wave conditions along the beach and in turns more uniform cross-shore processes. It is worth notifying that the temporal component of the translation mode is merely affected by the presence of the outcrop. However, the inclusion of the SLR effect leads to a reversal in the trend, which, combined to the shape of the spatial component, correspond to a progressive retreat of the coastline. While the SLR effect on the spatial component of the translation mode may be essentially related to the impact of variation of the mean sea level on wave transformation processes, the effect on the temporal component is likely to be the consequence of the application of the Bruun rule in LX-Shore.

Here, the assessment of the relevance of the new developments is based on an idealized test case without quantitative comparison to measurements. Future applications of the model will have to focus on assessing rigorously the model skill at real sites (e.g. [START_REF] Robinet | Controls of local geology and crossshore/longshore processes on embayed beach shoreline variability[END_REF] and the impact of calibration parameter uncertainties (e.g. D' Anna et al., in revision).

CONCLUSIONS

New developments were made in LX-Shore to further extend the range of possible real site applications. LX-Shore can now integrate real bathymetric features, which are shown to impact significantly the primary processes that drive shoreline variability. The model therefore appears as a promising tool to better understand and predict the respective contribution of the physical processes driving shoreline changes at sites exhibiting mixed (sandy/rocky) and complex seabed morphology. By extension, the model can also be used to anticipate the effects of building detached defense structures on shoreline variability, which are not necessarily confined to the vicinity of the structures. Finally, LX-Shore can now integrate SLR effects on the shoreline variability (through changes in wave transformation processes) and long-term retreat. This opens perspectives for shoreline prediction on the timescales from storm event to century.

  , black line); (ii) h r is kept offshore, i.e. for h i > D s (Figure 1, green line); (iii) a merged bathymetry (h m ) is computed and kept for D c ≤ h i ≤ D s to ensure a smooth transition, based on the formulas:

Figure 1 .

 1 Figure 1. Methodology used to merge the idealized bathymetry computed by LX-Shore with a measured bathymetry. Based on the idealized bathymetry, from 0 to Dc the idealized bathymetry is kept, from Dc to Ds the merging occurs and below Ds the measured bathymetry is kept. The bottom panel shows the sigmoid-shaped function used for the merging.

Figure 2 .Figure 3 .

 23 Figure 2. (a) Overview of the idealized embayed beach features and discretization adopted in LX-Shore. The red line denotes the initial planview shoreline used in all simulations. The dark grey patches denote the headland extents. The grid cell colors indicate the initial sediment

Figure 4 .

 4 Figure 4. Initial shoreline (dotted grey line) used in all simulations and 10-yr averages of the planview shorelines simulated in s1 (black line), s2 (blue line) and s3 (red line).

Figure 5 .

 5 Figure 5. Spatial (a,b,c) and temporal (d,e,f) components of the dominant modes of shoreline variability isolated by applying an EOF decomposition to the timestacks of planview shorelines simulated in simulations s1 (thick black line), s2 (blue line) and s3 (red line).

ACKNOWLEDGMENTS

The present work was funded by BRGM. The authors are grateful to the Water Research Laboratory team (University of New South Wales, Australia) which welcomed Arthur Robinet for a 3-month stay. Sydney wave data were kindly provided by the Manly Hydraulics Laboratory on behalf of the NSW Office of Environment and Heritage.

LITERATURE CITED

 Ashton, A.D. and Murray, A.B., 2006. High-


Abdel Halim

Melany Payoux

The link

come

The link between cognitive abilities and the main soft skills for the 21st Century workplace: A literature review

Abdel Halim., L & Payoux., M LPPL UPRES EA 4638, University of Nantes, France • Today's workforce is organized very differently than it was a few years ago.

Employees are required to adapt to changing environments, to communicate precisely (Sony & Mekoth, 2016), to think critically about tasks (Ozyurt, 2015) and to work in collaboration with others (Staar & Atzpodien, 2017).

• Therefore, to predict job performance, cognitive and soft skills often referred to as "21st century skills" are now central to recruitment, talent development and career management. • Also several studies have shown that we need at least seven cognitive abilities to thrive in the 21 st century workplace: working, semantic and episodic memory, verbal and numerical reasoning, selective attention and mental flexibility, verbal attention and inhibition, abstract reasoning, problem resolution and creativity and finally learning potential.

• Despite the acknowledged theoretical connection between cognitive skills and soft skills, such relationships have rarely been investigated and need to be examined because of its practical uses by HR professionals.

• So, the aim of this proposal is to examine the theoretical relationship between IDEAc's softs skills (Influence, Decision, Efficiency, Agility and Collaboration) and cognitive abilities previously mentioned.

1.Background

• We identified 22209 articles in our initial search.

• IDEAc's soft skills are covered with respectively those proportions: Decision (64,74%), Collaboration (27,44%), Influence (5,83%), Agility (1,97%) and Efficiency (0%).

2.Methods and Materials 4.Discussion & conclusions

• Preliminary results show that the link between IDEAc's softs skills and cognitive skills could be made.

• The lack of scientific studies which highlight the link between "Efficiency" and cognitive abilities (i.e., results) could be explained by the discrepancies between HR operational vocabulary on soft skills and those used in scientific studies. Redoing a literature review with other keywords could bring more significant results.

• It is, however, important to remain that we do not have a well-developed understanding of the causal direction, or the mechanisms underlying the relationships between soft skills and cognitive abilities.

• Nevertheless, this study can be useful in guiding future research design and analysis to these ends and what can also help HR professional to improve their decision-making skills.

• We conducted a broad literature search using the databases PsycINFO.

• We added one inclusion criteria : « 2009-2019 ».

• Search terms for IDEAc were « Influence », « Decision », « Efficiency », « Agility » and « Collaboration ».

• Search terms for the seven cognitive abilities were : «Working memory », « Semantic or Episodic memory », « Verbal and Numerical reasoning », « Attention », «Flexibility », « Abstract reasoning», and « Learning ».

• The search terms were combined in all 35 [5 (Soft skills IDEAc) x 7 (cognitive abilities)] possible ways.

• The review of the literature was limited to abstract.

3.Results

0

Contact Acknowledgment

We gratefully acknowledge the funding received towards our presentation from PerformanSe. 2016). The relationship between emotional intelligence, frontline employee adaptability, job satisfaction and job performance. Journal of Retailing and Consumer Services, 30, 20-32. • Trenvouez, A., Didry, A., Juhel, J., & Grasset, R. (2019). Soft skills for human resources professionals and behavioral styles in the workplace: A literature review. Poster presented at the WAPP 2019, 3rd World Conference on Personality, 2 -6, April, Hanoi, Vietnam.

•

  Özyurt, Ö. (2015). Examining the Critical Thinking Dispositions and the Problem Solving Skills of Computer Engineering Students. Eurasia Journal of Mathematics, Science & Technology Education, 11(2). • Staar, H., & Atzpodien, H. C. (2017). Vernetzte Zusammenarbeit erfolgreich gestalten-eine Analyse der Relevanz sozialer Fertigkeiten der Mitglieder in Gesundheitsnetzwerken. Gruppe. Interaktion.

Organisation. Zeitschrift für Angewandte Organisationspsychologie (GIO), 48(1), 25-40. • Sony, M., & Mekoth, N. (


H Michallet
email: herve.michallet@hmg.inpg.fr.

C Mathis

P Maïssa

F Dias
email: dias@cmla.ens-cachan.fr.

Flow Filling a Curved Pipe

A small scale experiment was designed to study the propagation of the front of a viscous fluid filling a curved pipe. Several Newtonian fluids with different viscosities and a non-Newtonian fluid have been used. The experiments show that there exists a minimum speed for completely filling the pipe, which depends on the parameters of the experiment (di-ameter d and radius of curvature R of the pipe, kinematic viscosity of the fluid). Appropriate dimensionless numbers are introduced to characterize the flow and optimal filling conditions.

Introduction

Prestressed concrete is commonly used in construction. To prevent the prestressed strands from any corrosion, a cement solution is generally injected in pipes of complex geometry. The pipe has to be completely filled in order to avoid weak zones in the structure. The aim of the present study is to observe experimentally the propagation of the front of a viscous fluid injected in a curved pipe, in order to understand the reasons for the appearance of unfilled zones. The solution used in construction is composed of cement, water, liquefier, and retarder. It is a non-Newtonian fluid with air bubbles and particles in suspension. As a first step, we mainly focus in this paper on the filling by Newtonian fluids.

Permanent flows in pipes have been widely studied ͓1͔, particularly in the case of curved pipes ͓2-4͔ and in engineering applications ͓5,6͔. However, pipe flows involving a propagating front or interface remain poorly understood. The problem of an airwater interface propagating in straight inclined pipes has been addressed only recently ͓7,8͔. As stated above, in industrial applications such as construction, the fluid is often mixed with air bubbles and the flow ought to be treated as a slug flow ͑see for example ͓9,10͔ for recent work on slug flows in pipes͒. But the emphasis of the present work is in the propagation of the interface which clearly separates the gas from the liquid.

The front propagation of a perfect fluid in a horizontal pipe was studied by Benjamin ͓11͔. See also the review by Simpson ͓12͔ on gravity currents for related studies and Asavenant and Vanden-Broeck ͓13͔. Benjamin analyzed the front in terms of a ''cavity flow'' displacing a fluid beneath it. The two-dimensional geometry is shown in Fig. 1͑a͒͑Benjamin also considered the case of a circular cross-section͒. Note that Benjamin studied only the case where the free surface detaches with an angle of 60 degrees. By applying conservation of mass, momentum, and energy, he found that there is a unique solution, characterized by

U ͱgd ϭ 1 2 , u ͱgh ϭͱ2, h d ϭ 1 2 ,
where the meaning of the various symbols is shown in Fig. 1. In other words, the cavity fills half of the box. An experiment was suggested by Benjamin to realize closely this flow. Liquid initially fills a long rectangular box closed at both ends and fixed horizontally. One end is then opened, and under the action of gravity the liquid flows out freely from this end. It can be expected that, after the transient effects of starting have disappeared, the air-filled cavity replacing the volume of the ejected fluid will progress steadily along the box. Observed in a frame of reference travelling with the front of the cavity, the motion of the liquid will appear to be steady, as shown in Fig. 1͑a͒. If the effects of viscosity and surface tension are neglected, the velocity of the cavity relative to a stationary observer will be U, and, since hϭd/2, the liquid will discharge from the open end with the same velocity. In addition to the solution studied by Benjamin, there is also a one-parameter family of solutions in which the free surface detaches from the box tangentially ͓14͔. Such a solution is plotted in Fig. 1͑b͒.A t the exit of the box, the Froude number u(gh) Ϫ1/2 for these solutions is between 1 and ͱ2. The limit ͱ2 corresponds to Ben- jamin's solution, while the limit 1 corresponds to the vanishing of the cavity. Using conservation of mass, momentum and energy, one finds the relation

u ͱgh ϭͱ d h . (1)
Unfortunately, simple arguments based on conservation of mass, momentum and energy do not give as much information in the case of a curved channel.

It is nevertheless important to keep in mind that the Froude number prevails to determine the size of the cavity and thus the ability to fill the pipe. This point will be discussed in the section where we consider various dimensionless numbers. Before that, we first describe the manner in which the experiments have been conducted. Then we present and discuss the experimental results.

Experimental Setup

The experimental setup is shown in Fig. 2. Gravity acts downwards. A transparent PVC pipe ͑''Tubclair''͒ of inner diameter d is curved with a radius R. Different diameters and radii were considered: 0.6рdр2c m ;5 рRр30 cm. The special case of a horizontal straight pipe was also investigated. The fluid was injected with a centrifugal pump at a constant mean flow velocity ͑noted U).

Several fluids were considered in order to investigate the effects of density , surface tension , and mainly kinematic viscosity ͑see Table 1͒. A series of runs has been performed with a non-Newtonian fluid ͑gel͒. Its viscosity has been measured with a plane-plane viscosimeter ͑see Fig. 3͒.

We assume that capillary effects are negligible, so that dӷL c where L c ϭ 1/2 (g) Ϫ1/2 is the capillary length (ϭ73 gs Ϫ2 for water and ϭ21 gs Ϫ2 for silicone͒. In the case of water L c ϭ0.26 cm and this assumption might be a shortcoming for small diameters. We nevertheless performed a few runs with water contaminated with wetting agents in order to diminish L c : the behavior and the results were the same within measurements errors. Moreover, the results presented in this paper are limited to dу1 cm.

Observations

For a very low flow rate ͑a typical mean flow velocity is U ϳ0.1 cm/s͒ and a sufficiently large pipe section ͑so that dӷL c ) the interface is more or less horizontal in the upflow part of the pipe. The flow is dominated by gravity. In the downflow part of the pipe the liquid is creeping down the pipe.

If U is slightly increased, the interface becomes more or less perpendicular to the pipe when rising. At the top, the interface deforms as shown in Fig. 2. The intersection between the interface and the upper part of the pipe ͑point C in Fig. 2͒ is a stagnation point.

If U is increased, the point C moves ͑interface in dotted line in Fig. 2͒. As soon as U is decreased back to zero, the pipe is emptying and the interface moves back to the top.

We consider that the pipe is completely filled when the point C has reached the point B ͑Fig. 2͒. We define U c as the critical speed beyond which the interface moves until the pipe is completely filled. The speed v is the mean velocity of the point C moving from A to B ͑Fig. 2͒:i fU is larger than U c , then v is greater than zero.

For the special case of the non-Newtonian fluid gel, the viscosity is smaller where the shear stress is larger, that is near the boundaries. Away from the boundaries, the behavior of the fluid is similar to a solid core pushed by the flow. In other words, the characteristic time of deformation of the interface is very long.

Measurements

The mean flow velocity U was obtained by measuring the flow rate in the steady-state regime, i.e., when the pipe was completely filled. The flow rate and v were estimated by measuring the filling time of a graded beaker and the travelling time of point C from A to B ͑Fig. 2͒ with a stopwatch. We may consider that these simple techniques do not lead to any bias error. By repeating the same run ͑at least five times͒ for different flow rates, we have estimated that the precision limit and thus the uncertainty was always better than 10 percent on U and v.

We measured v as a function of U for different experiments using various fluids, pipes, and curvatures ͑see Table 1͒. In Fig. 4, v is plotted versus U for several runs in four different configurations. The line vϭU delimits the region where physical solutions exist: v cannot be larger than U. In all our experiments v apparently varies linearly with U. We may thus determine the critical velocity U c ͑for vϭ0), the slope ␣ϭdv/dU and the interfacial velocity v(Uϭ0), by rms fitting a straight line into the data for each fluid, pipe diameter and radius of curvature.

Let us first consider the slope ␣ of the four straight lines in Fig. 4. The analysis presented in the Introduction considered a perfect fluid. It was assumed that the flow is identical in a reference frame moving with the front. In this framework, the filling of a pipe is similar to the emptying of the pipe. This would imply that ␣ is equal to one: adding ␦U to the mean entrance flow velocity would add the same value to the front propagation velocity. It is therefore not surprising that the slopes of the straight lines in Fig. 4 related to water are close to one. When the fluid is far from perfect but still Newtonian ͑case of the silicone oil V100 in Fig. 4͒, ␣ is clearly smaller than one: when the viscosity is large, the deformation of the front strongly depends on U. For a low viscosity, the deformation of the interface takes place on a very short time scale and the propagation of the front apparently does not depend on the reference frame. In the other limit, that is for a very large viscosity, the deformation of the front takes place on a very long time scale and the slope is again close to one ͑this is the case of the non-Newtonian fluid͒.

The velocity v(0) corresponds to the draining of the pipe. The measurements have been performed by bringing U to zero after filling the pipe: the flow is then from right to left in Fig. 2. It is important to emphasize that v(0) is very small in the case of the non-Newtonian fluid. We have estimated that the draining velocity is about 1 cm per week. The deformation of the front is visible only when the fluid is sheared and U c is therefore also close to zero for the gel. The agreement in v(0) between the measured values and the values deduced from the fitted straight lines in Fig. 4 is rather good except in the case of water for Rϭ30 cm. The discrepancy is probably due to measurement errors but also to wetting effects. The contact between the fluid and the PVC is indeed not the same when filling or emptying the pipe. It was not our purpose to investigate wetting effects and we did not test other materials for the pipe. We therefore have to consider that wetting effects may induce a bias error over U c of about 10 percent in the worst case, i.e., an uncertainty of 20 percent in the water case.

The critical speed U c ͑mean entrance velocity for which v ϭ0) is plotted against the pipe diameter d in Fig. 5. U c clearly tends to increase with d. Considering gravitational effects, the larger is d, the more efficient are gravity effects to deform the interface. It is thus not surprising that inertial effects have to be increased to compensate for it. Gravity effects could also reasonably depend on the radius of curvature. The interface deforms more rapidly when the pipe is curved: for a same diameter, U c is therefore smaller when Rϭϱ. Otherwise, if R is large but not infinite, the travelling time of the front in the pipe is large and the mean entrance velocity has to be large for the front to propagate. In other words, for a very small radius, the front has not enough time to deform before the exit of the pipe. This would explain the fact that U c is larger for Rϭ30 cm than for Rϭ5 cm.

In addition, above ϭ1 cm 2 /s, viscous effects tend to maintain the shape of the interface and thus decrease U c . This is noticeable for the silicone V100 and clearly visible for the silicone V300. Again the case of the gel requires specific comments. U c is very close to zero because the front deforms only when the fluid is sheared as mentioned above. The graphical estimation of U c gives a nonzero value only for dϭ1.5 cm. This configuration will allow us to estimate dimensionless numbers for the gel in the next section.

The data in Figs 4 and 5 are available in tabular form in Table 2.

Interpretation in Terms of Dimensionless Numbers

We may consider various dimensionless numbers that could characterize the flow in our special configuration. Note that in our experiments both the Bond number (Boϭd 2 g/) and the Weber number (WeϭdU 2 /) are very large: capillary effects are a lot weaker than inertial and gravitational effects.

Considering inertial and viscous effects in the framework of a curved pipe, we introduce a modified Reynolds number Re, which partly takes into account centrifugal effects, and is in fact the classical Reynolds number plus the Dean number ͓15͔:

Reϭ Ud ͩ1ϩͱ d R ͪ .
(2)

To balance inertial and gravity effects, we build a modified Froude number Fr, such as:

Frϭ U ͱgd ͩ1ϩͱ d R ͪ .
(3)

The modification of the classical dimensionless numbers may be understood by considering that there is a need to add centrifugal effects to maintain the fluid against the concave side of the pipe.

Note that ͑2͒ and ͑3͒ tend toward the classical Reynolds and Froude numbers in the limiting case of a straight pipe, i.e., when R tends toward infinity. Let us finally define K, which measures the balance between gravity and viscosity:

Kϭ

Fr Re ϭ ͱgd 3 .

(4) 1….

The straight lines are fitted curves by quadratic means.

Fig. 5

Critical speed versus the pipe diameter "see list of symbols in Table 1… Transactions of the ASME Note that K does not depend on the mean flow velocity U.Reand Fr are computed for UϭU c and plotted against K in Figs. 6 and7.

We may consider two kinds of estimation of the dimensionless numbers for the gel. First we may use the graphical estimation of U c from Fig. 4 for dϭ1.5 cm: U c ϳ1.2 cm/s, leading to ϳ200 cm 2 /s from Fig. 3 ͑assuming that ␥ ˙ϳ2U c /d). This allows us to plot the squares in Figs. 6 and7. Second, we may estimate from our measurement that U c ϳv(0) is nonzero but very small ͑six orders of magnitude less than the other values of velocity͒. In that case, the squares would be far off the bottom right-hand corners of Figs. 6 and7. In any case, the gel is possibly following the trend of Newtonian fluids with very large viscosity.

In Fig. 6, the data can be separated into two regions corresponding to partially ͑below the data͒ and completely ͑above the data͒ filled pipes. In a situation where the flow is characterized by a point belonging to the border line between the two regions, the filling would be effective over an infinite time, using an infinite amount of fluid. In the frame of a linear approximation, the data follow approximatively the (Ϫ1) slope, showing that Fr(U ϭU c ) is constant from Eq. ͑4͒.

Another way to represent the same results is given in Fig. 7. One sees that is of the same order of magnitude for all experiments while K varies over four orders of magnitude, but Fr is not rigorously constant. In the range of diameters we investigated, we can say roughly that Fr characterizes the ability to fill the pipes ͑Frϭ1 corresponds to the theoretical limit ͑1͒ where the flow detaches tangentially at the exit of a horizontal box, as suggested above in the introduction͒. Setting a flow rate leading to Frу1 would certainly fill any straight pipe.

A careful examination of Fig. 7 leads to several remarks. For Newtonian fluids, the dotted symbols, which correspond to straight pipes, lie below the points corresponding to curved pipes. The latter are gathered with no obvious dependence on the radius of curvature R. Therefore the modified Froude number, which takes into account centrifugal effects, can predict the effective filling of the curved pipes, but is of less interest in the limit of straight pipes. Each group of points corresponds to a different fluid, i.e., viscosity. The for Fr(UϭU c ) to increase with d is clearly visible in each group of points. For Kу10 Ϫ2 , Fr(U ϭU c ) is apparently decreasing with K and the flow would eventually tend to follow a non Newtonian behavior.

In almost all our experiments, Re is larger than 10 ͑see Fig. 6͒ indicating that viscous effects are smaller than inertial ones. Meanwhile, the data on Fr indicate that both inertial and gravitational effects are comparable.

Conclusion

In order to understand the failure in filling curved pipes, we performed a series of experiments using various fluids. We measured the velocity of propagation of the front v and we deduced from its linear dependence with the mean flow velocity U a critical value U c , below which there is no longer a complete filling of the pipe. This critical velocity U c is increasing with d for each fluid.

Results are presented in a synthetic way as relations between dimensionless numbers. We have shown, in the scope of our work, that a modified Froude number determines the ability of the flow to fill a curved pipe. In all our experiments and for Newtonian fluids of different viscosities, the filling was complete for Fr у0.6Ϯ0.2.

We also performed experiments with a non-Newtonian fluid whose rheological behavior is similar to the one of cement suspensions. The front deforms very slowly when the shear is weak because this increases the viscosity of the fluid. In that case, the complete filling may be obtained for Fr very small. Following ͓16,17͔, a more extensive study would be needed in order to understand the detail of the non-Newtonian effects.

For industrial applications we advise considering pipes of small diameters and injecting the solution at a relatively low flow rate in order to increase its effective viscosity. However, a future study of interest will be to repeat the present experiments in the context of slug flows.

Fig. 1

 1 Fig. 1 Steady two-dimensional flow past a cavity. "a… Solution with the free surface leaving the wall with a 60 degree angle. "b… Solution with the free surface leaving the wall tangentially. The point C denotes the detachment point.

Fig. 2

 2 Fig. 2 Diagram of the flow and notation. The interface is shown at two different times in solid and dotted lines.

Fig. 3

 3 Fig. 3 Dynamic viscosity of the gel versus the shear stress ␥ ˙. Its density is Ä1.03 gÕcm 3 . We note that this fluid is very viscous in the conditions of our experiments "we can estimate roughly that ␥ ˙È2UÕd is less than 10 s À1 …. Measurements have been performed with a plane-plane viscosimeter at Ecole des Mines de Paris "CEMEF, Sophia Antipolis….

Fig. 4

 4 Fig. 4 Speed of the front propagation v versus the mean entrance velocity U "dÄ1.5 cm; see list of symbols in Table1…. The straight lines are fitted curves by quadratic means.

Nomenclature d ϭ

 ϭ diameter of pipe R ϭ radius of curvature g ϭ acceleration of gravity ϭ kinematic viscosity ϭ density U ϭ mean entrance velocity u ϭ mean exit velocity v ϭ velocity of the upper point of the front U c ϭU(vϭ0) ϭ critical mean velocity h ϭ height of the exit flow FrϭU(1ϩ(d/R) 1/2 )/(gd) 1/2 ϭ modified Froude number ReϭUd(1ϩ(d/R)

Table 1 Fluids used in the experiments.

 1

		Density	Viscosity	Symbol
	Fluid ͑g/cm Water 1 Water ͑40 %͒ 1.16	0.01 0.086	᭹᭺ ᭢᭞	᭪ ᭞ •
	ϩ glycerol ͑60 %͒
	Silicone V50	1	0.5	᭤᭟
	Silicone V100	1	1	᭝
	Silicone V300	1	3	᭠
	Gel	1.03	͑S e eF i g .3 ͒	ᮀ

3 ) ͑cm 2 /s͒ Rϭ5c m Rϭ30 cm Rϭϱ

Table 2 Measured values

 2

									Table 2 "Continued…
			Radius	Pipe	Speed	Mean entrance			Radius	Pipe	Speed	Mean entrance
			of curvature	diameter	of front	velocity			of curvature	diameter	of front	velocity
		Fluid	R ͑cm͒	d ͑cm͒	v ͑cm/s͒	U ͑cm/s͒		Fluid	R ͑cm͒	d ͑cm͒	v ͑cm/s͒	U ͑cm/s͒
	᭤ Silicone V50 -	5 5	1.0 1.0	Ϫ8.6 0.6	0.0 10.2		--	ϱ ϱ	1.5 1.5	14.1 18.3	25.2 29.8
		-	5	1.0	3.2	14.1		-	ϱ	1.5	20.0	29.8
		-	5	1.0	11.8	28.7		-	ϱ	1.5	29.0	41.8
		-	5	1.0	18.9	36.6		-	ϱ	1.5	30.0	41.8
		-	5	1.2	Ϫ10.7	0.0	᭢ water ͑40%͒	5	1.2	1.6	16.7
		--	5 5	1.2 1.2	4.4 6.7	19.8 22.8	ϩ glyc. ͑60%͒ -	5	1.2	4.8	20.2
	᭟	--	5 30	1.2 1.0	12.2 Ϫ8.8	31.4 0.0	᭞	--	5 30	1.2 1.2	12.6 6.1	27.5 26.7
		-	30	1.0	1.0	13.0		-	30	1.2	15.0	38.8
		-	30	1.0	4.9	18.9		-	30	1.2	14.9	38.8
		-	30	1.0	8.9	23.6		-	30	1.2	29.5	63.2
		-	30	1.0	19.5	35.6		-	30	1.2	51.6	95.7
		-	30	1.0	22.1	38.4			30	2.0	5.7	39.3
		-	30	1.2	Ϫ10.8	0.0		-	30	2.0	3.3	42.4
		-	30	1.2	1.4	19.1		-	30	2.0	10.2	56.8
		-	30	1.2	4.5	24.0		-	30	2.0	15.2	64.2
		--	30 30	1.2 1.2	8.9 16.9	31.0 43.6	᭞ •	-	ϱ ϱ	2.0 2.0	0.0 0.0	6.9 8.2
	-----------⌬ silcone V100 ---------᭠ Silicone V300	30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30	1.5 1.5 1.5 1.5 1.5 2.0 2.0 2.0 2.0 2.0 2.0 1.2 1.2 1.2 1.2 1.2 1.5 1.5 1.5 1.5 1.5 1.5	Ϫ13.6 2.2 7.3 16.4 22.5 Ϫ18.0 2.0 3.4 4.8 6.4 11.2 Ϫ8.5 2.3 2.3 3.9 5.0 Ϫ11.0 1.3 6.8 10.3 11.8 Ϫ6.2	0.0 24.2 34.2 49.2 57.0 0.0 34.6 38.4 41.6 46.3 62.4 0.0 14.7 14.7 17.1 18.0 0.0 19.6 29.4 35.3 38.5 0.0	ᮀ gel	--------------------	ϱ ϱ ϱ ϱ ϱ 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30	2.0 2.0 2.0 2.0 2.0 1.2 1.2 1.2 1.2 1.5 1.5 1.5 1.5 1.5 1.5 2.0 2.0 2.0 2.0 2.0 2.0	16.4 22.7 34.9 45.3 57.0 0.0 1.8 20.3 22.3 0.0 6.2 6.5 11.1 13.7 42.3 0.0 38.7 35.9 31.0 44.1 55.9	24.8 32.1 47.1 69.2 80.4 0.0 2.0 23.1 27.8 0.0 9.3 9.6 15.1 18.5 46.0 0.0 49.6 53.6 55.6 64.0 88.9
		-	30	1.5	0.9	10.5
	᭹ water	5	1.0	Ϫ9.5	0.0
		-	5	1.0	3.0	12.4
		-	5	1.0	7.2	15.5
		-	5	1.0	10.8	18.1
		-	5	1.2	Ϫ10.3	0.0
		-	5	1.2	0.7	13.7
		-	5	1.2	1.4	14.0
		-	5	1.2	7.0	19.1
		-	5	1.5	Ϫ11.9	0.0
	᭺	-	30	1.0	Ϫ9.8	0.0
		-	30	1.0	3.2	14.8
		-	30	1.0	9.6	19.4
		-	30	1.0	13.5	23.1
		-	30	1.0	41.3	47.1
		-	30	1.2	Ϫ12.0	0.0
		-	0	1.2	0.5	16.9
		-	30	1.2	4.6	19.5
		-	30	1.2	5.6	22.5
		-	30	1.2	9.2	23.2
		-	30	1.2	14.1	29.0
		-	30	1.2	28.0	42.3
		-	30	1.2	28.8	44.9
		-	30	1.2	52.7	66.6
		-	30	1.5	Ϫ15.0	0.0
		-	30	1.5	9.9	28.9
		-	30	1.5	15.0	33.1
		-	30	1.5	14.9	33.1
		-	30	1.5	22.5	41.1
		-	30	1.5	32.6	49.2
		-	30	1.5	55.2	72.2
	ᶀ	--	ϱ ϱ	1.5 1.5	Ϫ12.0 0.0	0.0 9.4
		-	ϱ	1.5	0.0	13.3
		-	ϱ	1.5	7.1	17.8
		-	ϱ	1.5	9.2	20.6
		-	ϱ	1.5	13.0	25.2

computed for UÄU c versus K "4…" see list of symbols in Table 1.… Fig. 7 Modified Froude number "3… computed for U ÄU c ver- sus

  1/2 )/ ϭ modified Reynolds number K ϭ Fr/Re ϭ surface tension L c ϭ 1/2 (g) Ϫ1/2 ϭ capillary length Boϭd 2 g/ ϭ Bond number WeϭdU 2 / ϭ Weber number

	Fig. 6	Modified Reynolds number "2…

K "4…"see list of symbols in Table

1

.…

Acknowledgments

The authors thank C. Peiti at Ecole des Mines de Paris ͑CEMEF, Sophia Antipolis͒ for his help, and R. Leroy and O. Coussy for providing helpful comments. The authors are grateful to J.Ch. Bery for his help in the conception and realization of the experimental apparatus. This study was part of a contract with the Laboratoire Central des Ponts et Chausse ´es.


Stephany Mai Nishikawa

Soo Hyeon Kim

Yoshiho Ikeuchi

Timothée Levi

Teruo Fujii

NEURONAL CULTURE AT SINGLE CELL LEVEL FOR CREATION OF BIOLOGICAL NEURONAL NETWORK

Keywords: ITO electrode, Brain on a Chip, Cell control, Neuronal network

To understand the electrophysiological dynamics of neurons, the use of micro-electrode array (MEA) is common because it is non-invasive and it allows long-term monitoring. However, as the culture are randomly arranged on fixated electrodes, it is not capable of recording electrophysiology at single cell level. Here, we present a fully transparent neuronal culture device capable of monitoring cells growth or activity at single scale level. To trap cells singularly, dielectrophoresis is used on a ITO glass coated with fluoropolymer component called CYTOP® to manipulate cells.

INTRODUCTION

The use of idium tin oxide (ITO) has been employed for recording neuronal cell activity to investigate the biological phenomena such as dynamic neuronal function in a network. Microelectrodes array (MEA) is one of commercially available devices with ITO electrodes on glass that does not influence the healthy environment and it is usually applied in the investigation of neuronal activity experiment. Cell arrangement on such devices are usually randomized and cannot record each neurons. To understand the electrophysiological network pattern triggered by each neuron, there's a need to study single cell recording [START_REF] Yamamoto | A single-cell based hybrid neuronal network configured by integration of cell micropatterning and dynamic patch-clamp[END_REF]. Existing devices only records one by one or in small network of three. Here, a device capable of recording larger network of neurons with organized neuronal arrangement is concepted to understand cell dynamics in a patterned network, heterogeneity and single cell physiology. A fluoropolymer called CYTOP® which can make thin film coating of less than 1µm is employed for patterning and guiding the neuronal growth. Among the techniques of cell manipulation, dielectrophoresis (DEP) which uses electrokinetic phenomena is used for leading cell in a designated well. However, the use of DEP is very challenging for realizing neuronal network of in vitro culture due to high vulnerability to unfamiliar physiochemical environment such as material properties, humidity and pressure.

EXPERIMENTAL

Firstly, primary hippocampal neurons from mice embryo of D-18 were plated on a device with single channel made with CYTOP on ITO glass. After 10 days of maturation, neurons forming into an assembly were arranged in a single line with guided axonal projection (Fig1). Absence of neurons was also observed on CYTOP coated area. Secondly, a trap and drop system [START_REF] Kim | Efficient analysis of a small number of cancer cells at the single-cell level using an electroactive double-well array[END_REF] which was presented at MicroTas 2018 by Dr. Soo Hyeon Kim on microwellarray was reproduced for trapping the neurons singularly with DEP manipulation and dropping the cells to polydimethylsiloxane (PDMS) structured wells (Fig 2, Fig. 3). Larger wells made with PDMS material that rarely influence with neurophysiological environment were fabricated to enable spatial growth. As result, the cells were successfully trapped in the designated CYTOP wells of 10µm (Fig 4) and were dropped in the PDMS wells when the generator turned off.

RESULTS AND DISCUSSION

The use of thin coating of 1µm CYTOP successfully arranged the assembly of neurons in a single channel and guided the axonal projection on ITO glass electrode. Additionally, the completely transparent trap and drop device allowed single neuron trapping while having a full visibility though a microscope. However, only several cells were trapped and grew in the device and lacked axonal extension. To enhance the neuronal health, the device will be cultured with glial feeder layer in bilaminar model to support the neurons homeostasis and axonal projection [START_REF] Shimizu | Bilaminar co-culture of primary rat cortical neurons and glia[END_REF] The device was also innovated into a simpler design (without additional PDMS wells) with concepted channels between the wells to allow guided axonal projection (Fig 5).

CONCLUSION

Primary hippocampal neurons were manipulated successfully on a designated pattern on ITO electrode glass with 1µm patterned CYTOP. This single neuronal culture allowed cell monitoring for long term. An innovated neuronal network device with axonal projection channel is created and we hope to present the results on MicroTas 2019.

Figure 1 :

 1 Figure 1: Guided axon extension on patterned CYTOP Figure 2: Overview of the trap and drop system

Figure 3 :

 3 Figure 3: Concept of trap and drop device using dielectrophoresis technique. (A) primary neurons are attracted to the ITO glass coated with CYTOP patterned layer when electric field is applied. (B) primary neurons are dropped into a bigger well made with PDMS allowing (C) cells to grow.

Figure 4 .

 4 Figure 4. Primary neurons inside trap and drop device. Yellow arrows show single neural cells trapped in one well. Red arrow shows two different neuron cell sizes trapped in a well.

Figure 5 .

 5 Figure 5. (A) Transversal view of CYTOP wells on ITO glass. (B) Designed channels for 4 possible guided axonal extensions. (C) bi-directional axonal projection with larger well at the bottom of smaller cell trapping well for spatial growth


Stephany Mai Nishikawa
email: stephany@iis.u-tokyo.ac.jp

Farad Khoyratee

Soo Hyeon Kim

Yoshiho Ikeuchi

Kazuyuki Biomimetic Aihara

Teruo Fujii

Timothée Levi
email: levi@sat.t.u-tokyo.ac.jpwww.iis.u-tokyo.ac.jp

Spike

BIOMIMETIC SPIKE-TIMING BASED IONIC MICROSTIMULATION FOR NEURON CULTURE

Keywords: Ionic micro-stimulation, Neurons, Spiking Neural Network, Microfluidics

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Most of neurodegenerative diseases are incurable and debilitating conditions that influence cognitive and/or motor functions in millions of people worldwide. Neuroprosthesis are used today to support the quality of life but have yet to improve in their power consumption and biocompatibility issues. To realize a neuroprosthesis that are close to biological behavior, understanding the neurophysiological behaviors and the investigations on the interaction of neuronal cell assemblies is therefore essential. Often being silicon based, this kind of neuromorphic engineering is called biomimetic artificial neural systems [START_REF] Mahowald | A Silicon neuron[END_REF], [START_REF] Indiveri | Neuromorphic silicon neuron circuits[END_REF], [START_REF] Ambroise | Real-time biomimetic Central Pattern Generators into FPGA for hybrid experiments[END_REF]. The main goal of such systems is to design tools for biomedical applications including neuroprosthesis and to understand the human nervous system. This project reveals a new bio-hybrid system which includes a real-time Spiking Neural Network (SNN) and biomimetic ionic micro-stimulation coupled to living 'in vitro' neuron culture. To simplify the interactions with biological neurons and to reduce the bio-compatibility issues, we designed a biomimetic ionic micro-stimulation using microfluidic techniques and ionic exchange. Furthermore, to obtain a stimulation close to biological one, a neuromorphic system called Spiking Neural Network (SNN) is used. Firstly, the SNN was integrated in a digital platform FPGA. Then, the spike-timing coming from SNN system is used for triggering the microfluidic ionic micro-stimulation to make the stimulation the closest to biological behavior.

Biomimetic spike-timing based microfluidic system

There are three parts that compose the neuro-hybrid system (Fig. 1): (1) Spiking Neural Network (SNN), (2) neuron culture and (3) ionic micro-stimulation in microfluidic system. This system characterizes the neural network and its evolution by using biomimetic spiketiming based ionic stimulation. Each part will be described in detail throughout this article.

Fig. 1. Neuro-hybrid system with real-time digital SNN (green), biomimetic ionic micro-stimulation in microfluidic system (blue) and biological neuron culture (pink). SNN electrical output is converted into air pulse.

Spiking Neural Network (SNN)

The biomimetic SNN network is a neuromorphic system that is the closest detailed level of analogy to the nervous system [START_REF] Ambroise | Biorealistic spiking neural network on FPGA[END_REF]. A network of silicon neurons is connected via silicon synapses under plasticity rules. In neuro-hybrid experiments, spike timing and shape of the action potential (AP) must reproduce the same dynamics of a biological nerve impulse. This tunable biomimetic system works in real time and is based on Hodgkin-Huxley (HH) formalism [START_REF] Hodgkin | A quantitative description of membrane current and its applications to conduction and excitation in nerve[END_REF] which is the most bioplausible neuron model. We integrated the SNN in a digital platform FPGA. This spike timing-based system, SNN, was connected to the air pressure controller to trigger ionic micro-stimulation on the microfluidic device with stimulation closer to biological behavior.

Ionic micro-stimulation in microfluidic device

This device is structured with fluidic channels and pneumatic channels that allow the control of a stimulation sent to the neurons. The structure is described in Fig. 2 [START_REF] Nishikawa | Biomimetic microfluidic neurons for Biohybrid experiments[END_REF]. This device has four fluidic channels to input potassium chloride (KCl) of 10mM with each having a pneumatic valve to control the stimulation timing. The channels are equipped with ion selective permeable membrane (Sigma-Aldrich, Nafion® perfluorinated resin solution 5wt.%) to mimic the chemical exchange of biological neurons. It also has an integrated cell dish for the placement of living neurons with a thinness of less than 5mm allowing the observation through microscope during calcium imaging.

Neuron Culture

This bio-hybrid system contains 'in vitro' primary cultures of mice hippocampal neurons. Neurons were carefully dissected from embryos and were prepared as described [START_REF] Brewer | Optimized survival of hippocampal neurons in B27-supplemented Neurobasal, a new serumfree medium combination[END_REF]. Neurons isolated from embryonic day 17 ICR mouse were plated in Neurobasal medium. They were used in this experiment after two weeks of maturation inside the cell dish of the microfluidic device.

Results

Integration of Spiking Neural Network into FPGA

Regular and simplified equations of biomimetic SNN are implemented on a FPGA [START_REF] Levi | Digital implementation of Hodgkin-Huxley neuron model for neurological diseases studies[END_REF]. Electrical activities of different classes of cortical neurons can be simulated in biological time scale with strong correlation between hardware simulations (Fig. 3. red) and software simulations (Fig. 3. blue). The Central pattern generators (CPGs) provide bursts of spikes with frequency, variability and spike width in the ranges that are commonly observed in nature [START_REF] Ambroise | Biomimetic neural network for modifying biological dynamics during hybrid experiments[END_REF] Fig. 3. Comparison between model (blue) and FPGA (red). FPGA implementation and Matlab software simulations are nearly similar.

Micro stimulation in the microfluidic system

After two weeks of cell maturation, neuronal activities were visualized by calcium imaging with Fluo-4 AM (Molecular Probes, United States). To observe the spatial and temporal changes in calcium resulting from spontaneous activity, the raw sequences were processed to observe changes in fluorescence intensity of when the biomimetic ionic micro-stimulation is activated (Fig. 4).

Fig.

4.

Bio-hybrid experiments using biomimetic ionic micro-stimulation. When micro-stimulation is ON, the neuron culture is active. Calcium imaging technique is used for detecting neuronal activity.

When the 250ms length pulse of ON/OFF were sent to the pneumatic valve, the ionic stimulation was received by the neurons which in turn were stimulated by potassium ion. This was confirmed by staining the neurons with calcium imaging. The neural activities (Fig. 5.) successfully responded accordingly to the programmed ionic stimulation of 250ms. When the pneumatic microvalve is open, the fluorescent intensity increases. On the other hand, when it closes, the fluorescent decreases. During when the pneumatic valve is closed, the neurons activity returns back to its original intensity range.

Biomimetic spike-timing based

© The 2019 International Conference on Artificial Life and Robotics (ICAROB2019), Jan. 10-13, B-Con Plaza, Beppu, Oita, Japan Fig. 5. Stimulations were sent to the device with neuron culture. After calcium imaging, the activity was seen accordingly to the programmed air pulse sent to the pneumatic valve.

After this successful stimulation, the device was integrated with the Spiking Neural Network (SNN) to acquire a more bio-inspired stimulation. This was done by connecting the SNN implemented FPGA unit to the valve of the air controller.

Neural response to spike-timing based microstimulation

Fig. 6. After calcium imaging, the neural activity without any stimulation (top) and during the stimulation using SNN system (bottom) were recorded. The top graph shows the control activity of neurons at resting state. The bottom graph shows the neural activity when CPGs-based burst-like stimulations were sent for 500ms. Images were taken at 200ms frame rate and the fluorescent intensity mean was plotted. The neural activity is visualized when the florescent intensity increased spontaneously.

Activity was seen accordingly to the programmed spikes with added delays of approx. 250ms each spike.

The matured neurons were stained with calcium specific fluorescence to visualize the neural activity. The SNN stimulation sent has 8 bust-like spikes for 500ms with 500ms of no signals between stimulation. Fluorescence response observed shows a timely neural activity. However, there were delays in the response time for neuron activity (Fig. 6.). This may have been caused by the insufficient performance of the camera that was used.

Conclusion

Our device can generate a biomimetic stimulation accordingly to the sequenced air pressure programmed with Spiking Neural Network (SNN). The activity of neurons observed through calcium imaging showed timely response to SNN signals sent with certain delays.

The neural response was as expected when triggering the stimulation manually and by programmed air pressure of 250ms with air pressure controller. During the 500ms of OFF stimulation, the neural activities decreased as its fluorescent intensity diminished. On the other hand, as we integrated the Spiking Neural Network in the system, delays in response time for the activity of neuron after each stimulation were observed. However, during this recording, the computer had error for insufficient working space (RAM) and the recording showed some lags. This might have caused the error of response time in neural activities.

In the future, we will try to record the neural activity with a camera capable of capturing higher frame rate per second without having a lag during the recording of the neurons. On the other hand, the next step of our work is to create a closed-loop bio-hybrid system where the recorded neuronal activity is used as a feedback for SNN.

Fig. 2 .

 2 Fig.2. Structure of ionic micro-stimulation device. Set of four fluidic channels with under layer of pneumatic channel. The crossing of these two different types of channels will allow to control the input of solution in fluidic channel.


Damien Blanchard
email: damien.blanchard@u-bordeaux.fr

Kazuyuki Aihara
email: aihara@sat.t.u-tokyo.ac.jp

Timothée Levi
email: levi@sat.t.u-tokyo.ac.jpwww.iis.u-tokyo.ac.jp

Snake robot controlled by biomimetic CPGs

Keywords: snake robot, CPG, biomimetic, spiking neural network, locomotion

Locomotion is one of the most basic abilities in animals. Neurobiologists have established that locomotion results from the activity of half-center oscillators that provides alternation of bursts. Central Pattern Generators (CPGs) are neural networks capable of producing rhythmic patterned outputs without rhythmic sensory or central input. We propose a network of several biomimetic CPGs using biomimetic neuron model and synaptic plasticity. This network is implemented on a FPGA (Field Programmable Gate Array). We designed one unsupervised snake robot using this network of CPG. It is composed of one head wagon and 7 slave wagon. Infrared sensors are also embedded in the head wagon. This robot can reproduce the locomotion of one snake. The hardware robot is also described in this paper.

Introduction

Locomotion is one of the most basic abilities of animals. Neurobiologists have established that locomotion is created by Central Pattern Generator (CPG) activity that produces alternating bursts [START_REF] Brown | On the nature of the fundamental activity of the nervous centres; together with an analysis of the conditioning of rhythmic activity in progression and a theory of the evolution of function in the nervous system[END_REF] . Most rhythmic movements are programmed by CPG (Central Pattern Generator) networks [START_REF] Marder | Central pattern generators and the control of rhythmic movements[END_REF] . CPGs are neural networks capable of producing rhythmic patterned outputs without sensory input. CPGs can be found in animal locomotion such as swimming in salamander [START_REF] Ijspeert | From swimming to walking with a salamander robot driven by a spinal cord model[END_REF] , and lamprey [START_REF] Cohen | Modelling of intersegmental coordination in the lamprey central pattern generator for locomotion[END_REF] , as well as the cardiac pulsation system of the leech [START_REF] Gs | Bursting in Leech Heart Interneurons: Cell-Autonomous and Network-Based Mechanisms[END_REF][START_REF] Aa | A model of a segmental oscillator in the leech heartbeat neuronal network[END_REF] . Usually, in the field of robotics, CPGs are made using simple neuron models 7 or simple oscillators 8 , which are not biological time scales and which are considered as bio-inspired systems rather than biomimetic systems. These models provide sinusoidal oscillations and cannot reproduce biomimetic CPGs which have variability in performance. As our aim is to replicate biological behavior, we developed biomimetic CPGs using digital neuromorphic system.

Biomimetic CPGs in digital neuromorphic system

CPG neural network

The biomimetic CPG is based on the neural network system that controls the heartbeat of the leech [START_REF] Aa | A model of a segmental oscillator in the leech heartbeat neuronal network[END_REF] . This network is simple. It just only uses 8 excitatory neurons with 12 inhibitory synapses, making it an ideal candidate for our applications.

Simplification of the model

To design this CPG (Segmental Oscillator in Figure 1.C C), we use Hill's neural network diagram [START_REF] Aa | A model of a segmental oscillator in the leech heartbeat neuronal network[END_REF] . The network designed by Hill uses the Hodgkin-Huxley model [START_REF] Hodgkin | A quantitative description of membrane current and its applications to conduction and excitation in nerve[END_REF] and complex equations for synapses. A first work was to simplify this model and use Izhikevich neurons [START_REF] Izhikevich | Simple model of spiking neurons[END_REF] .

Using the short-term synaptic plasticity [START_REF] Izhikevich | Large-scale model of mammalian thalamocortical systems[END_REF] and the Izhikevich Regular Spiking neurons, we reproduced the activity of Hill's model.

Implementation in digital neuromorphic system

Biomimetic CPGs [START_REF] Barron-Zambrano | FPGA implementation of a configurable neuromorphic CPGbased locomotion controller[END_REF][START_REF] Li | Humanoids Learning to Walk : a Natural CPG-Actor-Critic Architecture[END_REF][START_REF] Donati | A spiking implementation of the lamprey's Central Pattern Generator in neuromorphic VLSI[END_REF][START_REF] Qi | A real-time FPGA implementation of a biologically inspired central pattern generator network[END_REF] are implemented in a digital platform: CMOD A7-A35 Xilinx Artix-7 FPGA (Field Programmable Gate Array). The neural network implementation architecture operates on a single computation core. This digital system works in real-time, requires few resources and is low power consumption [START_REF] Ambroise | Realtime biomimetic Central Pattern Generators in an FPGA for hybrid experiments[END_REF][START_REF] Ambroise | Biorealistic spiking neural network on FPGA[END_REF][START_REF] Grassia | Digital hardware implementation of a stochastic two-dimensional neuron model[END_REF][START_REF] Levi | Study of real-time biomimetic CPG on FPGA: behavior and evolution[END_REF] . Table 1 summarizes the resources. As we implemented the biomimetic CPGs in the FPGA board, we design the structure of the snake robot which embedded the FPGA board.

Snake robot

Design

The snake robots is composed of 7 wagons (like spinal cord) and 1 locomotive (mimics the brain). Biomimetic CPGs will control the different motors of the wagons. As the output of the FPGA board is 3.3V with low current, we need to add a power driver for the motors and one VHDL module for converting the CPGs to PWM (Pulse Width Modulation) to control motor speed. The PWM frequency should be between 100 Hz and 200 kHz.

The duty cycle of the PWM controls the speed of the motor. We use a Pololu_DRV8835 for the power driver.

Motors are DG01D with a 0.8 kg.cm couple and a no load speed of 90 rpm.

Fig. 4. 3D view of the wagon design

All the distance between the different wagons and the possible movement angle can be tuned thanks to mechanical structure. Figure 5 shows how the tuning is possible.

Fig. 5. The orange screw allows the tuning of distance of the different wagons but also the angle for the movement.

The design of electrical connection and components for the locomotive and for the wagon is described in figure 6.

Fig. 6. On the left, the electrical circuit of the wagon. On the right, the electrical circuit of the locomotive with the FPGA board.

The locomotive system has 3 IR sensors for detecting the different obstacles and then send stimulus to CPGs for modifying duty cycle and period. The tuning of the duty cycle allows the movement and the tuning of period allows the speed or the stop of the robot.

Results

Height CPGs were implemented (one for the locomotive and 7 for the wagons). Enable input are connected to each wagon to choose which ones are working for simulating a lesion of the spinal cord. The CPGs are connected for allowing forward locomotion like we can find in animal locomotion [START_REF] Matsuo | Biomimetic Motion Control System Based on a CPG for an Amphibious Multi-Link Mobile Robot[END_REF] . Fig. 7 describes this CPG network and the forward locomotion. The final snake robot is 67 cm long and works with realtime biomimetic CPGs and mimics the snake locomotion.

Figure 8 shows the final snake robot picture.

Fig. 8. Picture of the snake robot with the locomotive and the 7 wagons.

Conclusion

This article describes the design of a snake robot. This robot is controlled by real-time network of biomimetic CPGs. This system can be used in robotic applications [START_REF] Levi | Development and Applications of Biomimetic Neuronal Networks Toward BrainMorphic Artificial Intelligence[END_REF] for a closer behavior of biological animals. It can also be used in bio-hybrid robotics (biomimetic CPGs can be replaced by in vitro cell culture). We show previously that our system can be connected to biological neurons [START_REF] Joucla | Generation of locomotor-like activity in the isolated rat spinal cord using intraspinal electrical microstimulation driven by a digital neuromorphic CPG[END_REF][START_REF] Levi | Biomimetic neural network for modifying biological dynamics during hybrid experiments[END_REF] . Another advantage is that real-time systems can simulate multiple CPG topologies and simulate different hypothesis and protocol treatments for biomedical applications such as spinal cord injury [START_REF] Capogrosso | A brain-spine interface alleviating gait deficits after spinal cord injury in primates[END_REF] and neuroprosthetics [START_REF] Levi | Closed-loop systems for next-generation neuroprostheses[END_REF] .

Fig. 1 .

 1 Fig. 1. Electrical activity of the neural system of the Heartbeat of the Leech. Neurons are represented by circles. Inhibitory synapses are represented by black dots. A: Electrical Activity from the Leech [HILL 2003]. B: Diagram of an Elemental Oscillator of the Heart of the Leech. C: Diagram of the Segmental Oscillator of the Heart of the Leech.

Fig. 2 .

 2 Fig. 2. Comparison of CPG activity between Hill's model (in blue) and our model (in black) Even though we have simplified the model, we are able to reproduce the same behavior and characteristics as the biological CPGs. The important point is to allow a variability in the characteristics of the CPGs. Percentages of variability can be controlled by adjusting the parameters of short-term synaptic plasticity. By changing the parameters of the inhibitory synapses, it is possible to obtain different periods and cyclic ratios. This is very important especially for closed-loop experiments and to reproduce locomotion activities for robotics. A single parameter pair manages the modification of the period and of the duty cycle. The possible period interval in seconds is quite wide [0.3s; 35s].

Fig. 3 .

 3 Fig. 3. Three different CPGs with different period (1,7s, 3,6s and 5s). The output are from the scope.

Fig. 7 .

 7 Fig. 7. Forward locomotion with two examples of speed using a chain of 7 CPGs for the 7 wagons. Time axis is one second by division. Signal are the left neuron (N0, N2, N4, N6, N8, N10 and N12) of each CPG.

Table 1 .

 1 Utilization of the resources for the implementation of biomimetic CPGs

	Kintex-7	Utilization Available	Percentage
	LUT	7786	20800	37.43
	FF	3299	41600	7.93
	BRAM	50	50	100
	DSP	40	90	44.44
	Figure 3 describes the CMOD A7 FPGA board where
	CPGs with different period are implemented.

Acknowledgements

Authors thank Social Cooperation Program for Brain-Morphic AI to Resolve Social Issues, IIS, U Tokyo.


Yang Xia
email: xiayang@sat.t.u-tokyo.ac.jp

Seiji Uenohara
email: uenohara@sat.t.u-tokyo.ac.jp

Kazuyuki Aihara
email: aihara@sat.t.u-tokyo.ac.jp

Timothée Levi
email: levi@sat.t.u-tokyo.ac.jpwww.iis.u-tokyo.ac.jp

Timothée Levi Real-Time

Real-time implementation of ReSuMe learning in Spiking Neural Network

Keywords: Spiking neural network, ReSuMe, LIF, FPGA

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Neuromorphic systems are designed by mimicking or being inspired by the nervous system, which realizes * Typeset names in 10 pt Times Roman, uppercase. Use the footnote to indicate the present or permanent address of the author. robust, autonomous, and power-efficient information processing by highly parallel architecture. There are three common methods to realize the neuromorphic circuits, which are software 1 2 3 , analog hardware 4 5 6 7 and digital hardware 8 9 10 11 12 . Software can implement simple neuron model but a large scale neural network with complex neuron model cannot be realized in real-time. The power consumption is also quite important (kW for supercomputer).

For hardware implementation, compared to analog circuits, digital implementations consume more power but they are convenient to modify, more portable and lower cost for implementation with FPGA devices. Supervised learning was proposed as a successful concept of information processing in neural network [START_REF] Ghosh-Dastidar | A new supervised learning algorithm for multiple spiking neural networks with application in epilepsy and seizure detection[END_REF] . Recently, there has been an increasing body of evidence that instruction-based learning is also exploited by the brain. Remote Supervised Method (ReSuMe) is a new supervised learning method for Spiking Neural Networks.

The main reason for the study of ReSuMe is the need to invent an effective learning method to control the movement of people with physical disabilities. However, the in-depth analysis of ReSuMe method shows that this method is not only suitable for motion control tasks, but also suitable for other practical applications, including modeling, identification and control of various nonstationary and non-linear objects 14 15 . In this paper, we present a real-time ReSuMe learning implementation on FPGA using Leaky Integrate-and-fire (LIF) Spiking Neural Network (SNN). FPGA allows real-time implementation and embedded system [START_REF] Kechiche | Toward the Implementation of an ASIC-Like System on FPGA for Real-Time Video Processing with Power Reduction[END_REF] . We show that this implementation can make successful the learning on a specific pattern.

Method

This section proposed three methods that applied to the ReSuMe learning implementation on FPGA, which are LIF-neuron model, Postsynaptic potential (PSP) and Spike response model (SRM), as well as ReSuMe algorithm.

LIF neuron model

The LIF neuron is one of the simplest spiking neuron models. Due to the convenience with which it can be analyzed, simulated especially implemented in digital silicon neural network, the LIF neuron is very popular [START_REF] Gerstner | Spiking Neuron Models: Single Neurons, Populations, Plasticity[END_REF] . A neuron is modeled as a "leaky integrator" of its input I(t):

𝜏 𝑚 𝑑𝑣 𝑑𝑡 = -𝑣(𝑡) + 𝑅𝐼(𝑡) (1)
where v(t) represents the membrane potential at time t, 𝜏 𝑚 is the membrane time constant and R is the membrane resistance. This equation describes a simple resistorcapacitor (RC) circuit where the leakage term is due to the resistor and the integration of I(t) is due to the capacitor that is in parallel to the resistor. The spiking events are not explicitly modeled in the LIF model. Instead, when the membrane potential v(t) reaches a certain threshold 𝑣 𝑡ℎ (spiking threshold), it is instantaneously reset to a lower value 𝑣 𝑟 (reset potential) and the leaky integration process described by Eq. ( 1) starts a new with the initial value 𝑣 𝑟 . Consider the case of constant input: I(t) = I. We assume 𝑣 𝑟 = 0. The solution of Eq. ( 1) is then given by:

𝑣(𝑡) = 𝑅𝐼 [1 -exp (- 𝑡 𝜏 𝑚 )] (2)
Here v(t) is in an exponential decay. In discrete digital sequential circuit, a linear decay method is usually used to optimize computing process for saving hardware resources.

𝑑𝑣 = [-𝑣 + 𝑅𝐼] 𝑑𝑡 𝜏 𝑚 (3)
Eq. ( 3) describes the computing equation of dv, then solution v = v + dv obviously.

Postsynaptic Potential and Spike Response Model

By considering a single postsynaptic neuron i with a membrane potential 𝑢 𝑖 at time 𝑡 , a simplified SRM is defined [START_REF] Gerstner | Spiking Neuron Models: Single Neurons, Populations, Plasticity[END_REF] . 𝑢 𝑖 (𝑡|Χ, 𝑦 𝑖 ) ≔ ∑ 𝑤 𝑖𝑗 ∑ 𝜖(𝑡 -𝑡 𝑗 𝑓 )

𝑡 𝑗 𝑓 ∈𝑥 𝑗 𝑗 (4) This SRM signifies a dependence of the neuron's membrane potential on its presynaptic input pattern Χ from 𝑛 𝑖 synapses. An output spike occurs at a time 𝑡 𝑗 𝑓 .

The term of Eq. ( 4) describes a weighted summation of the pre-synaptic input: the 𝑤 𝑖𝑗 corresponds to the synaptic weight from a presynaptic neuron j, the kernel 𝜖 refers to the shape of an evoked PSP. The PSP kernel evolves according to

𝜖(𝑠) = 1 𝐶 ∫ exp (- 𝑠 ′ 𝜏 𝑚 )𝛼(𝑠 -𝑠′)𝑑𝑠′Θ(𝑠) ∞ 𝑆 ′ =0 (5)
The term Θ(𝑠)is the Heaviside step function defined such that Θ(𝑠) = 1 for 𝑠 ≥ 0 and Θ(𝑠) = 0. Here we approximate the postsynaptic current's time course by an exponential decay [START_REF] Gardner | [END_REF] .

𝛼(𝑠) = 𝑞 𝜏 𝑠 exp (- 𝑠 𝜏 𝑠 )Θ(𝑠) (6)
For a further simplified computer in digital circuit, the exponential decay 𝛼(𝑠) substitute the PSP kernel evolves 𝜖(𝑠) approximately.

ReSuMe architecture and algorithm

An implementation of ReSuMe in the Liquid State Machine (LSM) architecture is proposed as an example [START_REF] Maass | Computational models for generic cortical microcircuits[END_REF] .

The Liquid State Machine consists of a large, fixed "reservoir" network -the neural microcircuit (NMC) from which the desired output is obtained by training the suitable output connection weights.

In the implementation of ReSuMe method, the original LSM approach has been modified. The modified architecture consists of a set of input neurons 𝑁 𝑖𝑛 , the NMC structure, a set of learning neurons 𝑁 𝑙 with a total number k and a corresponding set of teacher neurons 𝑁 𝑑 (see Fig. 1). NMC receives signal 𝑠 𝑖𝑛 (𝑡) from 𝑁 𝑖𝑛 and transforms it into a vector of signals 𝑆 ̂𝑖𝑖𝑛 (𝑡) which i is presented to the learning neurons 𝑛 𝑖 𝑙 ∈𝑁 𝑙 . The teacher neurons 𝑁 𝑑 are not directly connected with any other structure.

Since we focus more on the ReSuMe learning implementation itself on this paper, so we generated presynapse to the learning neuron as NMC output.

The modification algorithm, which adjusts weights between pre-spike and post-neuron, is applied according to the following simplified equation:

The Fig. 2 shows a specific weights update process.

Implementation

These section proposed implementation of above methods we introduced above with results showed in waveforms.

Implementation of LIF neuron

We implement LIF Neuron with VHDL language in FPGA. By adjusting the size of dt and matching different time constant, our LIF Neuron can work at very high clock frequency (10 kHz), which means that its calculation accuracy is very high and the real-time requirement is realized.

Implementation of PSP and SRM

As we introduced in 2.2, we use simplified exponential decay to achieve RSM and PSP. Because post neurons are connected to 500 pre-synaptic inputs, hardware resources are still unacceptable if 500 exponential operations (even in linear decay) are performed in the same clock cycle. We adopt time division multiplexing, and use two-stage pipeline to complete 500 sets of PSP operations in 1000 clock cycles (actually 501 cycles, remaining standby, theoretically supporting input ceiling of 999), just using one multiplier and one adder.

Implementation of ReSuMe learning

ReSuMe we implemented includes an exponential attenuation (linear attenuation) which attribute to change of parameter k. Each time teacher input spike or postsynaptic neuron spike arrives, weight is updated. 500 exponential operations are performed. We use 500 clock cycles, time division multiplexing to achieve this change, and update weight in real-time so that PSP and RSM modules are used to calculate the correct weight instantly.

Architecture of ReSuMe learning

The overall hardware architecture of ReSuMe learning is shown in Fig. 6. We use LIF Neuron as post-synaptic neuron, equipped with ReSuMe learning module, and 500 pre-synaptic inputs are connected to post-synaptic neuron. Each connection is operated by PSP and summarized by RSM.

Results and discussion

Firstly, we simulate the actual learning waveform of ReSuMe. As shown in Fig. 7, we have completed the learning of 500 input data containing target patterns and random patterns. Post neurons can spike the location of each pattern with a minimal delay in 2-4 time steps. Here time step means cycle of computing and input patterns. For example, if clock frequency is 100 MHz and computing needs 1000 cycles under our implementation, the time step is 0.01 ms. Then we tested different input pattern groups, and completed real-time test and verification on the FPGA.

Conclusion

This paper introduces the advantages of using devices such as FPGA to realize the digital neural network. Methods of LIF neuron, PSP and RSM module as well as ReSuMe learning are described and illustrated on hardware implementation. Then the overall framework of ReSuMe learning is elaborated, with the different output data due to different inputs analyzed.

  𝑡) = [𝑆 𝑑 (𝑡) -𝑆 𝑑 (𝑡)] * [𝑎 𝑑 + ∫ 𝑊 𝑑 (𝑠 𝑑 )𝑆 𝑖𝑛 (𝑡 -𝑠 𝑑 )𝑑𝑠 𝑑 ∞ 0 ]

Fig. 2 .

 2 Fig. 2. An example waveform of synaptic weight updating in ReSuMe learning process.

Fig. 3 .

 3 Fig. 3. Stimulation of a LIF neuron by a constant input current: the time-course of the membrane potential, v(t) with 15 mV threshold and 20 mA current I (left); f-I curve for a LIF neuron (right).

Fig. 1 .

 1 Fig. 1. ReSuMe implemented in the modified Liquid State Machine Architecture.

Fig. 4 .

 4 Fig. 4. Waveform of PSP process in which indicates once presynaptic input spikes, the PSP starts an exponential decay.

Fig. 6 .

 6 Fig. 6. Architecture of ReSuMe learning.

Fig. 5 .

 5 Fig.5. Waveform of ReSuMe weight adjusting process which achieves the ideal functions as showed in Fig.2

Fig. 7 .

 7 Fig. 7. Waveform of different input patterns and their ReSuMe learning process. The third result shows that the minimum spike interval is 5 time steps.

© The 2019 International Conference on Artificial Life and Robotics (ICAROB2019), Jan. 10-13, B-Con Plaza, Beppu, Oita, Japan

Acknowledgements

Authors would like to thank Prof. Takashi Kohno, from IIS, the University of Tokyo.


Stephany Mai Nishikawa

Farad Khoyratee

Zhongyue Luo

Toshiharu Shiraishi

Kazuyuki Aihara

Yoshiho Ikeuchi

Soo Hyeon Kim

Teruo Fujii

Timothée Levi

NEURO-HYBRID SYSTEM WITH SPIKING NEURAL NETWORK AND BIOMIMETIC IONIC MICRO-STIMULATION

Keywords: Neuro-hybrid, Spiking Neural Network, micro-stimulation, microfluidic, biomimetic

Millions of people worldwide have neurodegenerative diseases that influence one's cognitive and/or motor functions. To bring neuroprosthesis into realization and for future long-term replacement of damaged brain areas with artificial devices, investigations on the interaction of neuronal cell assemblies are essential. To circumvent the limitations, we propose a new bio-hybrid system which includes a real-time Spiking Neural Network (SNN) and a biomimetic ionic micro-stimulation coupled to living in vitro neuron culture. This system characterizes the neural network and its evolution by using biomimetic spike-timing based ionic stimulation.

INTRODUCTION

Research on neurological disorder recently relies on alternative treatment using biomimetic and real-time hardware known as neuroprosthesis. Realization of such devices must follow several requirements including the knowledge on the relation between cells and machine or the neurons activities. This kind of neuromorphic engineering is called biomimetic artificial neural systems which are mostly in silicon-based [START_REF] Mahowald | A Silicon neuron[END_REF], [START_REF] Indiveri | Neuromorphic silicon neuron circuits[END_REF], [START_REF] Ambroise | Real-time biomimetic Central Pattern Generators into FPGA for hybrid experiments[END_REF]. The main goal of such systems is to design tools for biomedical applications including neuroprosthesis and to understand the human nervous system. This project reveals a new bio-hybrid system composed of a real-time Spiking Neural Network (SNN) and biomimetic ionic micro-stimulation coupled to living in vitro neuron culture. We describe here the three parts of this neuro-hybrid system (Figure 1). This system characterizes the neural network and its evolution by using biomimetic spike-timing based ionic stimulation.

Figure 1. Neuro-hybrid system with real-time digital SNN (green), biomimetic ionic micro-stimulation in microfluidic system (blue) and biological neuron culture (pink). SNN electrical output is converted into air pulse

THEORY

The biomimetic SNN network is a neuromorphic system that is the closest detailed level of analogy to the nervous system [START_REF] Ambroise | Biorealistic spiking neural network on FPGA[END_REF]. A network of silicon neurons is connected via silicon synapses under plasticity rules. In bio-hybrid experiments, spike timing and shape of the action potential (AP) must reproduce the same dynamics of a real nerve impulse. This tunable biomimetic system works in real time and is based on Hodgkin-Huxley (HH) formalism [START_REF] Hodgkin | A quantitative description of membrane current and its applications to conduction and excitation in nerve[END_REF] which is the most bio-plausible neuron model. We integrated the SNN in a digital platform (FPGA) (Figure 2). The spike timing of our SNN is used for triggering the microfluidic ionic micro-stimulation to make the stimulation the closest to biological behavior. Using silicon neurons bring some bio-compatibility issues (reject, power consumption). To simplify the interactions with biological neurons and to reduce the bio-compatibility issues keeping the same ionic currents (K+) to the biologic ones, we designed a biomimetic ionic micro-stimulation using microfluidic techniques and ionic exchange. This microfluidic device [START_REF] Levi | Microfluidic Neurons, a New Way in Neuromorphic Engineering?[END_REF] is made with Polydimethylsiloxane (PDMS) which is highly biocompatible and it is composed of two chambers intra-and extracellular modelling connected by channels with different KCl concentrations actuated by pneumatic microvalves [START_REF] Unger | Monolithic Microfabricated Valves and Pumps by Multilayer Soft Lithography[END_REF] (Figure 3). The pneumatic valve-based technology requires multilayer soft lithography and an external air pressure controller. The latter air pressure system is triggered by the SNN to obtain biological like neuronal behavior. The channels are equipped with ion selective permeable membrane to mimic the chemical exchange of biological neurons.

EXPERIMENTAL

This bio-hybrid system contains 'in vitro' primary cultures of mice hippocampal neurons. Neurons were carefully dissected from embryos and were prepared as described [START_REF] Brewer | Optimized survival of hippocampal neurons in B27-supplemented Neurobasal, a new serum-free medium combination[END_REF]. Neurons isolated from embryonic day 17 ICR mouse were plated in Neurobasal medium. The experiments were proceeded after two weeks of maturation inside the cell dish of the microfluidic device. Neuronal activities were visualized by calcium imaging with Fluo-4 AM. To visualize the spatial and temporal changes in calcium resulting from spontaneous activity, the raw sequences were processed to highlight changes in fluorescence intensity of when the biomimetic ionic micro-stimulation is activated (Figure 4).

RESULTS

Regular and simplified equations of biomimetic SNN are implemented on a FPGA [START_REF] Levi | Digital implementation of Hodgkin-Huxley neuron model for neurological diseases studies[END_REF]. Electrical activities of different classes of cortical neurons can be simulated in biological time scale. Mathematical tools show strong correlation between hardware simulations and software simulations (Figure 5). The device with integrated cell dish for the placement of living neurons has a thinness of less than 5mm allowing the observation through microscope during calcium imaging. Our device can generate a biomimetic stimulation [START_REF] Nishikawa | Biomimetic microfluidic neurons for Bio-hybrid experiments[END_REF] accordingly to the sequenced air pressure. The opening of the pneumatic microvalve was sequentially repeated after three pulses of 0.1s, 0.25s and 0.5s with each separated with 0.2s of closure. The matured neurons were stained with calcium specific fluorescence to observe the neural activity. Fluorescence response was successfully observed only during the opening of pneumatic microvalve.

CONCLUSION

We demonstrated successful micro-stimulations of the neuron cells with potassium ions using microfluidic pneumatic microvalves operated by air pressure controller. The neural activities responded according to the programmed air pressure applied (0.1s, 0.25s and 0.5s), capable of generating a biomimetic behavior. Next step is to perform bio-hybrid experiments coupling all the modules. This system will investigate neural network dynamics and test new neuro-prosthesis [START_REF] Levi | Closed-Loop Systems for Next-Generation Neuroprostheses[END_REF], [START_REF] Hochberg | Neuronal ensemble control of prosthetic devices by a human with tetraplegia[END_REF], [START_REF] Nicolelis | Principles of neural ensemble physiology underlying the operation of brainmachine interfaces[END_REF], [START_REF] Bonifazi | In vitro large-scale experimental and theoretical studies for the realization of bi-directional brain-prostheses[END_REF], [START_REF] Ambroise | Biomimetic neural network for modifying biological dynamics during hybrid experiments[END_REF].

Figure 2 .

 2 Figure 2. Real-time SNN. HH neuron model are implemented in a FPGA. Regular spiking neuron model parameters are: neuron surface area of 96 µm², a K, Na and leakage conductance respectively at 5.00, 50.00, and 0.07 mS.cm² and a K, Na and leakage reversal potential respectively at 50, 100, -70 mV. The SNN works in biological time and reproduces the spike timing of biological neural network.

Figure 3 .

 3 Figure 3. Different layers of ionic micro-stimulation device. In this example, four micro-stimulation holes are designed.

Figure 4 .

 4 Figure 4. Bio-hybrid experiments using biomimetic ionic micro-stimulation. When micro-stimulation is ON, the neuron culture is active. Calcium imaging technique is used for detecting neuronal activity.

Figure 5 .

 5 Figure 5. Real-time hardware implementation on FPGA of different neuron family. On the left panel, inhibitory neurons named Fast Spiking neuron. On the right panel, excitatory neurons named Regular Spiking neuron. FPGA implementation and Matlab software simulations are similar.

Figure 6 .

 6 Figure 6. Calcium Imaging of opening (ON) and closing (OFF) of pneumatic microvalve

ACKNOWLEDGEMENTS

This work was partially supported by NEC Corporation


Stefano Buccelli

Ilaria Colombi

Yannick Bornat

Valentina Pasquale

Laura Martines

Francesco Difato

Alverna Alberto

Jacopo Tessadori

M Tedesco

Paolo Massobrio

Yannick

Implementing a hybrid bi-directional communication between neuronal assemblies

come   L'archive ouverte pluridisciplinaire

I. INTRODUCTION

The consequences of brain lesions represent a major economic and healthcare priority worldwide. One of the main challenges is to enhance the recovery of motor and cognitive functions following such insults. Sometimes pharmacological treatments and physical rehabilitation are not enough to regain a fully independent life [START_REF] Soekadar | Brain-machine interfaces in neurorehabilitation of stroke[END_REF]. Recent technological improvements in the field of neural interfaces are opening new frontiers to meet the growing demand for alternative solutions [START_REF] Guggenmos | Restoration of function after brain damage using a neural prosthesis[END_REF]. In vitro studies can provide an easy and low-risk model to explore promising intervention strategies.

Here we realized an in vitro platform for next generation brain-prostheses aimed at artificially link two previously connected neuronal assemblies.

II. METHODS

The experimental set-up consisted of a biological element and an artificial system (i.e. the neuroprosthetic device). We Tecnologia, Genova, 16163, Italy (e-mail: Ilaria.colombi@iit.it); Y. B. is with IMS Lab, University of Bordeaux, Bordeaux, 33400, France (e-mail: yannick.bornat@ims-bordeaux.fr); V. P. is with NBT, Istituto Italiano di Tecnologia, Genova, 16163, Italy (e-mail: valentina.pasquale@iit.it); L. M. is with the NBT, Istituto Italiano di Tecnologia, Genova, 16163, Italy (e-mail: martines.laura@gmail.com); F. D. is with NBT, Istituto Italiano di Tecnologia, Genova, 16163, Italy (e-mail: francesco.difato@iit.it); A. A. is with with Rehab Technologies IIT-INAIL joint Lab, Istituto Italiano di Tecnologia, Genova, 16163, Italy (e-mail: alberto.averna@iit.it); J. T. is with ADVR, Istituto Italiano di Tecnologia, Genova, 16163, Italy (e-mail: jacopo.tessadori@iit.it); M. T. is with DIBRIS, University of Genova, 16145, (email: brunella.tedesco@unige.it); P. M. is with DIBRIS, University of Genova, 16145, Italy (email: paolo.massobrio@unige.it); T. L. is with IMS Lab, University of Bordeaux, 33400, France (e-mail: timothee.levi@u-bordeaux.fr); M. C. is with with Rehab Technologies IIT-INAIL joint Lab, Istituto Italiano di Tecnologia, Genova, 16163, Italty (corresponding author, phone: +39 010 71781; fax: +39 010 71781 230; e-mail: michela.chiappalone@iit.it). used dissociated cortical rat neurons plated over a 60-channel Micro-Electrode Array (MEA) characterized by two interconnected networks (i.e. bimodular), exhibiting spiking and network bursting (NB) activity. During experiments, we performed a laser ablation to cut the axonal connections between the two neuronal modules to mimic a focal lesion.

The neuroprosthetic device consisted of a neuromorphic board based on Field Programmable Gate Array (FPGA) and was able to perform real-time spike/NB detection. Upon detection, an activity dependent stimulation paradigm was applied from one module to the other and vice versa (Bi-directional Bridging, BB). When an artificial Spiking Neural Network was configured in the FPGA, we were able to implement a Hybrid Bidirectional Bridging (HBB) experimental paradigm, where the SNN was replacing one of the two modules of the biological network.

III. RESULTS

We first characterized the effect of a focal lesion on a bimodular network. After a few hours of recording, the firing of the network reached pre-lesion values thus proving the viability of the network. The correlation between the modules collapsed to zero proving both an anatomical and functional disconnection. We then implemented two neuroprosthetic strategies: BB and HBB. Through the BB strategy, we were able to partially recover the cross correlation features without imposing any preferred direction. The HBB protocol was performed when we aimed at replacing one entire module by means of a SNN. We used NB events to control the stimulation and also in this configuration, we obtained a partial recovery of the cross correlation without affecting the firing of the survived module.

IV. CONCLUSIONS

Our neuroprosthetic device is able to artificially restore the bi-directional interaction between two previously connected neural ssemblies.. The developed stimulation strategies can be seen as a starting point towards a new rehabilitation paradigm that can act directly at the CNS level.

  Italy (e-mail: stefano.buccelli@iit.it); I. C. is with with Rehab Technologies IIT-INAIL joint Lab, Istituto Italiano di

	*Research supported by European Commission, ICT-FET
	FP7/2007-2013, FET Young Explorers, grant agreement N° 284772 BRAIN
	BOW (www.brainbowproject.eu).
	S. B. is with Rehab Technologies IIT-INAIL joint Lab, Istituto Italiano di
	Tecnologia, Genova, 16163,