Mehmet Kévin Pons
email: kevin.pons@principia.fr

Mehmet Ersoy
email: mehmet.ersoy@univ-tln.fr

Adaptive mesh refinement method. Part 1: Automatic thresholding based on a distribution function

Adaptive mesh refinement method. Part 1: Automatic thresholding based on a distribution function

Introduction

A major issue of many modelling challenges is to solve accurately processes over very large ranges in spatial scale. Solving accurately those processes with high resolution inexorably leads to heavy computational time [START_REF] Popinet | Quadtree-adaptive tsunami modelling[END_REF][START_REF] Popinet | Adaptive modelling of long-distance wave propagation and fine-scale flooding during the tohoku tsunami[END_REF][START_REF] Altazin | Numerical investigation of bb-amr scheme using entropy production as refinement criterion[END_REF].

In principle, Adaptive Mesh Refinement (AMR) [START_REF] Berger | Adaptive mesh refinement for hyperbolic partial differential equations[END_REF] method allow to solve in a reasonable CPU time equations adjusting the computational effort locally to maintain a quite uniform level of accuracy. It relies on macro cells which can be refined (and then possibly coarsened).

The zones where the mesh needs to be refined are determined thanks to a suitable mesh refinement criterion and a mesh refinement threshold. Depending on the application, equations and numerical methods, a variety of criteria might be used based on an error estimation procedure or a feature detection technique. Basically, the zones where a given threshold is exceeded are refined. Usually, the mesh refinement threshold is a test case dependent parameter. If this parameter is set too small, the results will be accurate but the computational cost will be expensive. On the contrary, if this threshold is set too large, pertinent regions are not detected and the mesh adaptation algorithm will not be efficient. Generally speaking, if this parameter is not well-calibrated, it may yield to unnecessary refined cells, and sometimes the numerical solution can develop spurious oscillations. Therefore, the way to choose the mesh refinement threshold is very important.

In this paper, we propose a method to set automatically the threshold parameter to detect almost all relevant regions to refine. The construction of the automatic selection of the threshold is based on the decreasing rearrangement function of the mesh refinement criterion. The decreasing rearrangement provides a description of the mesh refinement criterion in terms of the local maxima which are sorted from the smallest to the largest. Thus, in principle, this method allows to set the threshold automatically small enough to detect the most relevant cells to refine. In practice, we are able to reach a good balance between the accuracy and the computational time without any hand-calibration of the mesh refinement threshold as usually performed in mesh adaptation algorithms.

In the sequel, we focus on h-AMR method. As a non-exhaustive list of references, one can refer to [START_REF] Berger | Adaptive mesh refinement for hyperbolic partial differential equations[END_REF][START_REF] Berger | Local adaptive mesh refinement for shock hydrodynamics[END_REF][START_REF] Verfürth | A posteriori error estimation and adaptive mesh-refinement techniques[END_REF][START_REF] Garcia | Adaptive mesh and algorithm refinement using direct simulation monte carlo[END_REF][START_REF] Kirk | libmesh: a c++ library for parallel adaptive mesh refinement/coarsening simulations[END_REF][START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF].

The paper is organised as follows. The Section 2 is devoted to the presentation of the finite volume method and the principle of the h-AMR approach. We also present some mesh refinement criteria used in this paper. In Sect. 3, we discuss about advantages and drawbacks of some classical thresholding methods. Then, we introduce a new method based on the decreasing rearrangement of the mesh refinement criterion to construct a suitable automatic threshold. Finally, for a given one dimensional Riemann problem for the non-linear shallow water equations, we illustrate the efficiency of the method to detect simple waves without hand-calibration. In the second part of this paper [START_REF] Pons | Adaptive mesh refinement method. Part 2: Application to tsunamis propagation[END_REF], more complex multi-dimensional test cases are considered.

Finite volume approximation and h-AMR algorithm

For the sake of simplicity, we consider here the one-dimensional case (we refer to the second part of this paper [START_REF] Pons | Adaptive mesh refinement method. Part 2: Application to tsunamis propagation[END_REF] for the multi-dimensional case and applications). This section summarises the main features of the finite volume approximation and the adaptive mesh refinement method for a general non-linear hyperbolic system

∂ t w + ∂ x f(w) = 0 w(x, 0) = w 0 (x) (1)
where x ∈ R is the coordinate, t > 0 is the time, w(x,t) ∈ R d is the unknown vector with d ≥ 1, f(w) ∈ R d is the flux at a space-time point (x,t) and w 0 is the initial data.

In particular, in this work, we will consider the Saint-Venant equations as a prototype model with w = h hu and f(w

) = hu hu 2 + g h 2 2
. Here the unknowns h(x,t) ≥ 0 and u(x,t) ∈ R stand respectively for the water height and the depthaveraged velocity of the water at a space-time point (x,t) and g ≈ 9.81 m/s 2 is the gravitational constant. To illustrate the method of automatic thresholding proposed in Sect. 3.2, we will consider the Riemann problem (c.f. Sect. 3.3):

w 0 (x) = w l if x < 0 w r if x > 0 (2)
Interested readers can found multi-dimensional test cases in the second part of this paper [START_REF] Pons | Adaptive mesh refinement method. Part 2: Application to tsunamis propagation[END_REF].

First order finite volume approximation

We suppose that the computational domain

Ω ⊂ R is split into a set of cells C k = x k-1/2 , x k+1/2 of size δ x k with x k±1/2 = x k ± δ x k /2 such that Ω = ∪ k C k .
We assume that k ∈ Z to simplify the presentation. We define the discrete time by t n+1 = t n + δt n where the time step δt n satisfies a Courant, Friedrichs, Levy (CFL) condition. On a given cell C k of center x k , noting

w w w k (t) 1 δ x k C k w w w (x,t) dx .
the approximation of the mean value of the unknown w w w(x,t) on C k at time t n , integrating (1) over C k × (t n ,t n+1 ), and dividing by δ x k , we obtain the first order finite volume approximation of Eqs. (1) (see for instance [START_REF] Godlewski | Numerical approximation of hyperbolic systems of conservation laws[END_REF][START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics[END_REF][START_REF] Eymard | Finite volume methods[END_REF]):

w n+1 k = w n k - δt n δ x k F F F k+1/2 (w n k , w n k+1 ) -F F F k-1/2 (w n k-1 , w n k )
where

F F F k+1/2 (w n k , w n k+1 ) is the approximation of the flux F F F k+1/2 (w n k , w n k+1 ) ≈ 1 δt n t n+1 t n f(w(x k+1/2 ,t)) dt . F F F k+1/2 (w n k , w n k+1
) is defined via any three-point scheme (see for instance [START_REF] Perthame | A kinetic scheme for the Saint-Venant system with a source term[END_REF][START_REF] Simeoni | Remarks on the consistency of upwind source at interface schemes on nonuniform grids[END_REF][START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics[END_REF]). In this paper, we consider a first order space-time discretisation using the Godunov scheme, i.e. the numerical flux F F F k+1/2 is computed with the exact solution of the 1D Riemann problem at the interface x k+1/2 with the left state w n k and the right state w n k+1 (for further details see, for instance, [START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics[END_REF]).

Principle of h-AMR method

For the convenience of the reader, we recall the procedure for the first order scheme in the one dimensional case. We refer to [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF] for the second order scheme and to [START_REF] Altazin | Numerical investigation of bb-amr scheme using entropy production as refinement criterion[END_REF][START_REF] Pons | Adaptive mesh refinement method. Part 2: Application to tsunamis propagation[END_REF] for the multi-dimensional case. Let us note that in h-AMR method the mesh evolves in time, thus each cell C k should be written C n k . To simplify the notation, we omit here the time dependency.

In order to reduce the time necessary to manage the refinement, we use "macro cells" which could be refined by generating hierarchical grids. Each cell can be split in two sub-cells. We thus produce a dyadic cells graph, whose numbering (in basis 2) allows a quick computing scan to determine the adjacent cells. We use the following notations: let k b be the index associated to the macro cell numbered k with b a binary number which contains the hierarchical information of a sub-cell. In particular, the level of a sub-cell C k b is defined as the length(b) -1. By convention, the coarsest cells are of level l = 1, while the finest cells are of level l max ∈ N * . For instance, a macro cell C k 0 of level 1 will be split into two sub-cells C k 00 and C k 01 of level 2. A mesh refinement example is proposed in Fig. 1. The mesh refinement procedure can be simply expressed as follows. We introduce a non negative mesh refinement criterion S n k b computed on each cell C k b at time t n and compared, for instance, to the average

S m S m = 1 |Ω | ∑ k b S n k b (3)
This approach is the so-called mean method (c.f. Sect. 3.1 for more details and a presentation of some other approaches).

We then define two coefficients 0 < β min ≤ β max , which determine the ratio of the cells to be refined or coarsened. Thus, for each cell C k b :

• if S n k b > α max = S m β max ,
the cell is refined and split into two sub-cells C k b0 and C k b1 , • if S n k b0 < α min = S m β min and S n k b1 < α min , the cell is coarsened into a cell C k b . Therefore, the adaptive algorithm stops when the level l of a cell C k b reaches the maximum level l max .

Remark 1. The threshold parameters β min and β max allow to set a percentage of mesh refinement and mesh coarsening with respect to the quantity S m . The grid is not refined or coarsened between α min and α max . It is not surprising that these settings will deteriorate or improve the accuracy of the numerical solution. For instance, the more β min and β max is small, the more accurate are the results to the expense of the computational time. In many AMR method, the choice of these thresholds parameters is therefore a critical weakness.

On one hand, if a cell C k b is split into two sub-cells C k b0 and C k b1 , averaged values at time t n are projected on each sub-cell:

w w w n k b0 = w w w n k b1 = w w w n k b
On the other hand, if two sub-cells C k b0 and C k b1 are coarsened, we initialise the new cell C k b as:

w w w n k b = 1 2 w w w n k b0 + w w w n k b1
Moreover, as done in [8, Section 4.2], we use a numerical smoothing grid technique which prevents two adjacent cells to have a level difference greater than two for stability purpose.

To simplify the notations, we do not use the binary subscript notation and we will assume that β min = β max = β and therefore α min = α max = α for the sake of simplicity. By abuse of notations, we write the adaptive algorithm as follows: for each cell C k :

• if S n k > α, the cell is refined and split into two sub-cells, • if S n k < α, the cell is coarsened (in the sense that two cells C k b0 and C k b1 are coarsened into a cell C k b as described before).

If the threshold is not correctly chosen, it is well-known that such a strategy may yields to spurious oscillations since locally, the grid incessantly changes between coarse and fine cells. However, even in this settings (i.e. α min = α max = α), we will see that the method proposed for the automatic selection of the threshold does not develop, in general, spurious oscillations. If the threshold parameter is automatically well-chosen at each time step, one can always define thereafter α min and α max to allow to set a percentage of mesh refinement and mesh coarsening to reinforce the stability of the algorithm.

Mesh refinement criteria

The choice of a suitable mesh refinement criterion is crucial to perform accurate mesh adaptation. Indeed, ideally, one would like to set up a mesh that minimises the discretisation error and maximises the computational efficiency of the algorithm. The construction of an efficient criterion remains a difficult subject, since it depends on physical phenomena and on numerical methods. Multi-resolution analysis [START_REF] Harten | Multiresolution algorithms for the numerical solution of hyperbolic conservation laws[END_REF][START_REF] Müller | Adaptive multiscale schemes for conservation laws[END_REF] provides some robust tools to estimate and control rigorously the error at each refinement level. Based upon Richardson-type estimates of the local truncation error [START_REF] Berger | Adaptive mesh refinement for hyperbolic partial differential equations[END_REF], one can also refine or coarsen to reach a given accuracy for a minimum amount of work. However, near discontinuous solutions, Richardson extrapolation is invalid [START_REF] Powell | Adaptive-mesh algorithms for computational fluid dynamics[END_REF]. It will work but it is no more reliable than heuristic methods and the extra cost cannot be justified. The entropy indicator for hyperbolic systems is also an efficient mesh refinement criterion for smooth and non smooth solutions [START_REF] Puppo | Numerical entropy and adaptivity for finite volume schemes[END_REF]. Finally, the mesh criterion can be also constructed, as simple as possible, as the combination of heuristic gradients based on physical quantities to detect pertinent regions to refine or to coarsen [START_REF] Quirk | A parallel adaptive grid algorithm for computational shock hydrodynamics[END_REF]. We also refer to other methods based on a posteriori error estimates based on a measure of the local residual [START_REF] Verfürth | A posteriori error estimation and adaptive mesh-refinement techniques[END_REF] and [START_REF] Karni | A smoothness indicator for adaptive algorithms for hyperbolic systems[END_REF].

In this paper, we will consider the following three mesh refinement criteria.

Criterion 0: Since the exact solution of the Riemann problem (1)-( 2) is known, we consider the criterion as the "best 1 mesh refinement criterion"

S n k := |h n k -h ex (x k ,t n )| (4)
where h ex stands for the exact solution. Criterion 1: One can show that System (1) admits a mathematical entropy E E(w) = q 2 2h + gh 2 2 with q = hu which satisfies the entropy inequality

E (x,t) = ∂ t E(w(x,t)) + ∂ x (G(w(x,t))) ≤ 0 . (5)
with

G(w) = E(w) + gh 2 2 u .
The entropy satisfies a conservation equation only in regions where the solution is smooth and an inequality when the solution develops shocks. In simple cases, it can be proved that the term missing in [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF] to make it an equality is a Dirac mass.

Following [START_REF] Puppo | Numerical entropy and adaptivity for finite volume schemes[END_REF][START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF][START_REF] Altazin | Numerical investigation of bb-amr scheme using entropy production as refinement criterion[END_REF], we define the numerical density of entropy production 1 The criterion 0 is qualified as "the best" in the a posteriori sense. Indeed it gives the exact numerical error at time t n . However it will not necessary give the best convergence rate at the end because the best mesh refinement should be done where large errors are going to be done (at time t n+1 ).

E n k = E(w n+1 k ) -E(w n k ) δt n + G k+1/2 (w n k , w n k+1 ) -G k-1/2 (w n k-1 , w n k ) δ x k
where the numerical flux G k+1/2 is computed with the exact solution of the 1D Riemann problem at the interface x k+1/2 with the left state w n k and the right state w n k+1 as done before in Sect. 2.1. At the discrete level, as shown in [START_REF] Puppo | Numerical entropy and adaptivity for finite volume schemes[END_REF], E n k converges to zero with the same rate of the local truncation error in smooth regions. Therefore, E n k can be regarded as an a posteriori error indicator. Moreover, as proved in [START_REF] Puppo | Numerical entropy and adaptivity for finite volume schemes[END_REF], the entropy production built on first order monotone schemes is essentially negative definite, in the sense that positive overshoots in the cell entropy inequality are possible, but only on non monotone data and their amplitudes decrease fast under grid refinement for smooth solutions. Thus, in the following, we define

S n k = |E n k | . (6)
Criterion 2: Finally, the last mesh refinement criterion is based on the gradient method. It can be constructed as a combination of heuristic gradients based on physical quantities to detect pertinent regions to refine or to coarsen [START_REF] Quirk | A parallel adaptive grid algorithm for computational shock hydrodynamics[END_REF]. To make it as simple as possible, we consider the gradient of h which can be an efficient mesh refinement criterion if the threshold is well-chosen,

S n k = h n k+1 -h n k x k+1 -x k . (7)
In what follows, we define the piecewise mesh refinement criterion as follows :

S(x,t) = ∑ k S n k 1 C k ×[t n ,t n+1 [ (x) (8)
where

1 A (x) = 1 if x ∈ A , 0 otherwise.

Mesh refinement thresholding methods

Generally speaking, as already described before, the mesh adaptation process is constructed as follows: if S(x,t) > α, the cell is refined and split into two sub-cells and otherwise the cell is coarsened. The quantity α is the so-called mesh refinement threshold. As presented before, if we set α = β S m for some β > 0, we obtain the mean method. We propose here to review some of the classical thresholding methods. Then, we present a new automatic thresholding approach constructed with the decreasing rearrangement function of the mesh refinement criterion S.

A brief review of the classical thresholding methods

In general, the mesh refinement criterion of error indicator type can be written as a function of the local discretisation error ε

S(x,t) = f (ε(x,t))
where f : R → R + . If the mesh refinement criterion is an a priori or a posteriori error estimates, then the function f is in general linear as for the criterion (4), say S(x,t) = κε(x,t), with κ > 0. In this settings, it is easy to define the mesh refinement threshold α = ε(x,t)/κ which almost corresponds to the desired local accuracy.

Let us note that in general it is very difficult to construct such a mesh refinement criterion. The criterion ( 6) is an a posteriori error indicator. In the case of the mesh refinement criterion ( 6), the function is not necessarily linear and therefore it is not easy to define a suitable threshold. There exists several thresholding methods. We recall the classical ones: mean method, mean method with standard deviation, two-steps filtering, wavelets, local maxima, Dannenhoffer and Powell cumulative distribution methods and we finally present our method. The last one is the only method, up to our knowledge, yielding to an efficient automatic threshold almost parameterless contrary to the previous ones.

Mean method.

It is the simplest and the cheapest method. It can be an efficient method to detect homogeneous smooth solutions, i.e. solutions for which the local maxima of S have almost the same amplitudes. In what follows, we refer to smooth solutions for homogeneous smooth solutions. As done before, for a given time t > 0, the mesh refinement criterion S is compared to

α = β 1 |Ω | Ω S(x,t) dx
where β is a user-calibrated dimensionless parameter. For a discontinuous solution, the mean method is not able to detect the relevant smooth regions to refine if the threshold β is not set small enough. However, the mean method is, in general, suitable for smooth solutions for which setting β = 1 yields to an efficient refinement. In Fig. 2, we have illustrated formally a discontinuous solution, see Fig. 2(a), and a smooth one, see Fig. 2(b). In Fig. 2(a), the peak localised around at x = x B = 3.75 shows the presence of a discontinuity of the solution while the remainder region corresponds to the smooth part of the solution.

On one hand, for a discontinuous solution, one can see in Fig. 2(a) that if the parameter β ≥ 1 (i.e. α ≥ S m ) then the smooth part of the solution will not be marked for refining. On the other hand, for a smooth solution, one can observe in Fig. 2(b) that setting β = 1 allows to detect relevant part of the mesh to refine. In general, to be efficient β should be calibrated such that β ∈ (0, 1] or equivalently α ∈ (0, S m ]. Mean method with standard deviation.

The mean method can be improved by taking into account the fluctuations of the mesh refinement criterion. For a given time t > 0, the mesh refinement criterion S is now compared to

α = β 1 |Ω | Ω S(x,t) dx + δ σ (x,t)
where σ is the standard deviation of S(x,t) and δ ∈ R. These methods are largely used due to their simplicities, see for instance [START_REF] Kallinderis | Adaptation methods for a new navier-stokes algorithm[END_REF][START_REF] Warren | Grid convergence for adaptive methods[END_REF][START_REF] Jasak | Automatic resolution control for the finite-volume method, part 2: Adaptive mesh refinement and coarsening[END_REF][START_REF] Schmidt | Influence of adaptive mesh refinement and the hydro solver on shear-induced mass stripping in a minor-merger scenario[END_REF] however β and δ are problem dependent and therefore user-calibrated parameters.

Filtering: two steps method.

The method of "filtering" is constructed to detect the smooth part of a discontinuous solution, see for instance [START_REF] Warren | Grid convergence for adaptive methods[END_REF][START_REF] Aftosmis | Upwind method for simulation of viscous flow on adaptively refined meshes[END_REF]. This method improves the classical mean method by applying two successive filtering. For a given time t > 0 and for a given fixed β , as a first step, the mean method is applied to S to identify discontinuities and steep gradients in the solution, like in the example in Fig. 2(a). Then, we redefine the mesh refinement criterion as

S 1 (x,t) = S(x,t) if S(x,t) < α = β 1 |Ω | Ω S(x,t) dx, otherwise S 1 (x,t) = β 1 |Ω | Ω S(x,t) dx.
We now apply the mean method to the mesh refinement criterion S 1 . As a consequence, this method allows to filter the large peaks into the function S and try to detect also the smooth part of the solution. However, the two steps may be not enough to efficiently capture the smooth part of the solution. Moreover, the second drawback concerns still the parameter β which is again test-case dependent.

Filtering: wavelet method.

In contrast with the previous method, the detection of smooth and discontinuous solutions can be efficiently captured with a space-time localisation based on a wavelet transform, see for instance [START_REF] Leonard | A wavelet-based adaptive mesh refinement criterion for large-eddy simulation[END_REF]. A nice feature of wavelets are the detection of discontinuities and high transitions along with the resulting absolute values of wavelet coefficients are large. The localisation of such a region is generally well-captured for smooth regions2 as displayed in Fig. 3(b) for a given signal

x ∈ [0, 1], S(x) = exp(-1000(x -0.7) 2 ) + 3 exp(-5000(x -0.2) 2 ) + 100 if x ∈ [0.4, 0.401] 0.3r 0,1 (x) if x ∈ [0.401, 0.43]
where x → r 0,1 (x) ∈ [0, 1] returns a random number, see Fig. 3(a). High variations of the signal yield to large absolute values of wavelet coefficients centered around the discontinuity at all regions. Depending on the wavelet's support, the larger the region is, the larger the set of coefficients affected in the wavelet transform is. It defines the so-called "cone of influence". As a consequence, the discontinuity has the smallest region. On the contrary, smooth signal produces relatively large wavelet coefficients at large regions and again the definition of the cone of influence holds. The interpretation of the absolute value of wavelet coefficients as a mesh refinement criterion is based on the cone of influence and its support. For instance, in the example illustrated in Fig. 3, for a given region (which corresponds to a threshold parameter α), say y = 100, the intersection with the cone of influence provides three intervals located around the point x = 0.2, x = 0.4 and x = 0.7 which are the regions to refine. However, as for the method described before, we are still confronted with the choice of a threshold parameter (y). Moreover, even if this method is based on fast wavelet transform, it increases the global cost of the numerical method. Local maxima approach.

It consists in looking for the local maxima of the mesh refinement criterion. As a first drawback, the local maxima of the mesh refinement criterion are not necessarily those of the local error. Once a local maximum is detected, say x * , the surrounding area to refine should be large enough. It can be defined without parameter by constructing the interval

[x * -ζ (x * ), x * + ζ (x * )] where ζ (x *
) is the distance between the local maxima position x * and the closest inflexion point of the criterion. As a second drawback, such a method can be quite costly in a multi-dimensional settings. Moreover, for a noisy criterion function (like for instance in Fig. 5(a)), it can lead to an over refined domain.

In Fig. 4, we illustrate the method for a given mesh refinement criterion. The intervals to refine around each local maxima are defined by the use of the closest inflexion point (see Fig. 4(b) and Fig. 4(a)). Dannenhoffer cumulative distribution method and Powell method.

The method proposed by Dannenhoffer [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF] is simple and a priori efficient to automatically set the threshold parameter.

Let us consider a given discrete mesh refinement criterion (8) for a given time t > 0. Let (α j ) 0≤ j≤M be an increasing sequence of M + 1 threshold parameters such that α 0 = 0 and α M = max x S(x,t n ) where α j , is given by

α j = S M j M γ , 0 ≤ j ≤ M (9)
with γ ≥ 1. Following Dannenhoffer [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF], we construct the piecewise constant distribution function d j = d(α j ) as:

d(α) = M-1 ∑ j=0 d j 1 (α j ,α j+1 ) (α) with d j = #{k ; S n k > α j } ( 10
)
where # is the number of elements in the set {S n k > α j }. This function is useful to detect the local maxima.

We display the shape of the cumulative distribution function in Fig. 5(b) (black area) for a given perturbed mesh refinement criterion (see Fig. 5(a)) as in [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF]. In this illustration, the first and the second peaks around x = 0.5 and x = 0.7 can be associated to a steep gradient or a discontinuity in the solution. The peak around x = 0.25 can be associated to a region where the solution is smooth. Ideally, we would like to detect those cells to mark for refinement. One can define several threshold parameters denoted α = A, α = B, α = C and α = D in Fig. 5. The threshold α = A corresponds to the area for which the mesh refinement criterion oscillates. Thus, setting α = A yields to detect almost all the cells to refine and it is not a good strategy. The threshold D is too large and almost all regions are not detected. Thus, Dannenhoffer proposes to set a threshold between α = C and α = B which is localised at the "knee" of the cumulative distribution (see Fig. 5(b)). In particular, following [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF], the "knee" of the cumulative distribution is defined as the smallest α = α Dan such that

d (α) = -1 .
Ideally, if the location of the "knee" can be computed then one can construct an automatic threshold which allows to refine the pertinent areas. However, as mentioned in [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF], to compute correctly the threshold α Dan , a construction of an accurate interpolation of the derivative is required which increases the computational cost of the AMR algorithm. Otherwise, the threshold α Dan computed is often too small and the method detects too many cells to refine as mentioned in [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF].

Based on this cumulative distribution function, Powell [START_REF] Powell | An embedded mesh procedure for leading-edge vortex flows[END_REF] proposes to automatically set the threshold at the lowest value of the refinement parameter that produces a local maximum in the curvature κ = |d | (1+d 2 ) 3/2 . The threshold α = α P is localised before the first inflexion point of the cumulative distribution function but still an accurate interpolation of the second and third derivatives are required.

As a conclusion, in general, the construction of a suitable mesh refinement threshold is either user-calibrated or automatically but with an extra numerical cost. In the following section, we propose a new automatic mesh refinement threshold α PE , only based on the cumulative distribution function and not on its derivatives. The threshold α PE is defined such that max (α Dan , α P ) < α PE . The computational cost of the construction of the threshold α PE is less expensive than the Dannenhoffer [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF] or the Powell [START_REF] Powell | An embedded mesh procedure for leading-edge vortex flows[END_REF] method.

A distribution function for automatic thresholding

In this section, we present the new method to select the threshold automatically to detect relevant cells to refine.

For the sake of simplicity and readability, we present our new threshold method in the one-dimensional case and we refer to the second part of this paper [START_REF] Pons | Adaptive mesh refinement method. Part 2: Application to tsunamis propagation[END_REF] for the generalisation to the multi-dimensional case.

Let us consider a smooth (at least twice differentiable) mesh refinement criterion S(x,t) ∈ R + , x ∈ [0, L] and t > 0 where L is the length of the domain. The time t being fixed, we write in the sequel S = S(x). Without loss of generality, we suppose that S(0) = S(L) = S (0) = S (L) = 0 and 0 < S ∞ = max x∈(0,L) S(x) < ∞ (if S ∞ = 0 the numerical solution is identically equal to zero).

Remark 2. The meaning of the assumption S(0) = S(L) = S (0) = S (L) = 0 depends on the mesh refinement criterion and simplifies the following computations. For the mesh refinement criterion [START_REF] Ersoy | A simple and efficient new algorithm to increase the regularity and vanishing moments of biorthogonal wavelets[END_REF], it means that the free surface is flat at the boundaries. Generally speaking, if for instance at x = 0, S(0) = 0 and S (0) = 0, one can always arbitrarily extend the function S on the interval (-ε, 0) such that S(-ε) = S (-ε) = 0 for a given ε > 0. So without loss of generality, one can consider that S(0) = S(L) = S (0) = S (L) = 0.

In view of the above assumptions, the set

Z α = {x ∈ (0, L); ϕ α (x) = S(x) -α = 0 and S (x) = 0}
is not empty. Indeed, for each 0 < α < S ∞ , since S ∞ > 0, S has at least one maximum. Then, there exists p α ∈ N * such that the number of elements in the set Z α is

#Z α = 2p α .
Thus, for all α ∈ (0, S ∞ ), one can describe the set Z α as follows

Z α = {x 0 (α) < x 1 (α) < • • • < x 2p α -2 (α) < x 2p α -1 (α)} .
Let us assume that S has p local maxima. Then there exists an increasing sequence α * k 1≤k≤p and a sequence x * k 1≤k≤p such that

∀k = 1, . . . , p S (x * k ) = 0, S (x * k ) < 0 and S(x * k ) = α * k .
By definition, the sequence α * k 1≤k≤p represents the local maxima from the smallest to the largest with α * p = S ∞ . Let us also remark that

x * k / ∈ Z α , ∀k = 1, . . . , p, ∀α ∈ (0, S ∞ ).
With these settings, we define the cumulative distribution function

α ∈ [0, S ∞ ] → d(α) :=    L if α = 0 , ∑ p α k=1 x 2k+1 (α) -x 2k (α) if 0 < α < S ∞ , 0 if α = S ∞ . ( 11
)
Remark 3. The function [START_REF] Godlewski | Numerical approximation of hyperbolic systems of conservation laws[END_REF] corresponds to the well-known decreasing rearrangement function of S which is widely used in optimal transport problems. For α ∈ [0, S ∞ ], d(α) is the Lebesgue measure of the set {S(x) > α}.

Remark 4. From a theoretical viewpoint, in general, the presence of a Dirac mass in the mesh refinement criterion indicates a discontinuous solution. The automatic thresholding being based on the Lebesgue measure, the measure of a singleton being zero, the method will be not able to detect the presence of discontinuities. However, from a numerical viewpoint, due to the numerical diffusion, the Dirac mass is often associated to a steep peak and globally the mesh refinement criterion can be assumed smooth enough. As a consequence, the method works to detect large peaks which are associated to numerical approximation of discontinuous solutions. Therefore, we can assume that the mesh refinement criterion is at least twice differentiable.

Even if the smoothness assumption on S is rough, the method developed works to detect the discontinuity, the steep gradient and the smooth part of the solution.

Then, one has:

Theorem 1. Let l ∈ N * , S ∈ C l ([0, L], R + ) be a smooth function such that S(0) = S(L) = S (0) = S (L) = 0, 0 < S ∞ < ∞ with p local maxima. Then 1. S ∞ 0 d(α) dα = L 0 S(x) dx. 2. d ∈ C 0 ([0, S ∞ ], R + ) and a. ∀α ∈ [0, S ∞ ], d (α) < 0 b. ∀k ∈ 0, p , lim α→α * k d (α) = -∞ with the convention α * 0 := 0 3. d ∈ C l (D * , R + ) on the set D * := p-1 k=0 (α * k , α * k+1 ) .
Proof. By construction of the function [START_REF] Godlewski | Numerical approximation of hyperbolic systems of conservation laws[END_REF], the first and the second properties hold. For α ∈ [0, S ∞ ], the derivative of d is given by

d (α) = p α ∑ k=1 d dα x 2k+1 (α) - d dα x 2k (α) .
Since ∀k ∈ 0, p α , x 2k (α) ∈ Z α , i.e., ϕ α (x 2k (α)) = S(x 2k (α))α = 0, we get

d dα x 2k (α) = 1 S (x 2k (α))
.

Since S (x 2k (α)) > 0 and S (x 2k+1 (α)) < 0 by construction, we therefore deduce that

d (α) = p α ∑ k=1 S (x 2k (α)) -S (x 2k+1 (α)) S (x 2k+1 (α))S (x 2k (α)) < 0 .
The last property is obtained as a straightforward consequence of the above inequality.

In order to illustrate the previous notations and results, we display in Fig. 6, for a given mesh refinement criterion S (see Fig. 6(a)), the corresponding distribution function [START_REF] Godlewski | Numerical approximation of hyperbolic systems of conservation laws[END_REF] (see Fig. 6(b)), the first (see Fig. 6(c)) and the second derivative (see Fig. 6(d)). For better readability the distribution has been normed by L. In this example, we consider a function S with p = 3. In Fig. 6(b), the points for which the distribution function has vertical asymptotes are

α * 1 = S(x * 1 = 5) = 5, α * 2 = S(x * 2 = 2.5) = 10 and α * 3 = S(x * 1 = 7.5) = 20 = S ∞ .
and correspond to the local maxima of S which are sorted, by construction, from the smallest to the largest. In Fig. 6(a), for α = 2, the set Z α is composed of 6 = 2p α elements which are approximately

x 0 (α) ≈ 2.0957, x 1 (α) ≈ 2.9043, x 2 (α) ≈ 4.566,

x 3 (α) ≈ 5.434, x 4 (α) ≈ 7.3921 and x 5 (α) ≈ 7.6079 .

The distribution function at α = 2 is then computed as follows

d(α) = (x 1 (α) -x 0 (α)) + (x 3 (α) -x 2 (α)) + (x 5 (α) -x 4 (α)) L .
Fig. 6 The distribution function d, the first d and the second derivative d for the mesh refinement criterion

S(x) = L exp(-10(x -L/4) 2 ) + L/2 exp(-5(x -L/2) 2 ) + 2L exp(-200(x -3L/4) 2 ) with L = 10.
Remark 5. The cumulative distribution function [START_REF] Godlewski | Numerical approximation of hyperbolic systems of conservation laws[END_REF] is the continuous version of the Dannenhoffer cumulative distribution function [START_REF] Dannenhoffer | Grid adaptation for complex two-dimensional transonic flows[END_REF] (see also Eq. ( 10)).

In view of Theorem 1, the distribution function ( 11) is useful because it provides a complete description of the mesh refinement criterion S, in particular, it allows to localise the smallest and the largest local maximum. Therefore, the construction of an automatic threshold based on the function d can be relevant.

However, on one hand, the method proposed by Dannenhoffer (resp. by Powell) is based on the first derivative (resp. the second and the third derivatives) of the distribution function [START_REF] Garcia | Adaptive mesh and algorithm refinement using direct simulation monte carlo[END_REF] which requires an accurate interpolation. On the other hand, using only the distribution function [START_REF] Godlewski | Numerical approximation of hyperbolic systems of conservation laws[END_REF] or equivalently [START_REF] Garcia | Adaptive mesh and algorithm refinement using direct simulation monte carlo[END_REF] to construct a relevant mesh refinement threshold is numerically difficult since a local maximum corresponds to detect an area for which the derivative of d is infinite. Therefore, one way to overpass these drawbacks is to use the distribution function [START_REF] Godlewski | Numerical approximation of hyperbolic systems of conservation laws[END_REF], without using any derivatives of d, and simply weighting the function d as follows

f (α) = αd(α) . (12)
This function f transforms a local point of d between the first inflexion point of d and the first local maximum of S (i.e. the first local vertical asymptote of d) into a local maximum. More precisely, we have the following useful properties:

Corollary 1. Assume that S is twice differentiable and has p local maxima. Then, by virtue of Theorem 1, there exists

∀k = 0 . . . p -1, α * * k+1 ∈ (α * k , α * k+1 ) such that d (α * * k+1 ) = 0
and the function f (12) has p local maxima α 1 , . . . , α p such that

∀k = 1 . . . p, α k ∈ (α * * k , α * k ) .
Proof. The function f inherits the same regularity property of d as given in Theorem 1. Therefore if S is twice differentiable, then it is sufficient to show that on each interval (α * k , α * k+1 ), the function f is concave. The proof is similar for all k. Let us show the result for k = 1. For all α ∈ (α * * 1 , α * 1 ) since d (α) < 0 and f (α) = αd (α) + 2d (α) < 0, we deduce that f is a strictly concave function. As a consequence, there exists α 1 ∈ (α * * 1 , α * 1 ) which maximises f .

For each k, since α k ∈ (α * * k , α * k ), the local maximum α k of the function f does not coincide with the local maximum α * k of S. By construction, setting the threshold parameter α = α k for a given k yields to max (α Dan , α P ) < α k .

On one hand, as mentioned in Sect. 3.1, for smooth solutions, the mesh refinement threshold α = S m is in general relevant (S m = 1 |Ω | Ω S(x) dx). On the other hand, we have also pointed out that α should be less than S m for discontinuous solutions. Thus, the threshold α = α k if α k < S m , otherwise α = S m is a good candidate. More precisely, we define our automatic threshold as follows:

α PE such that f (α PE ) = max 0<α≤S m f (α) . (13)
One important consequence of this choice is that the threshold is a global maximum of a regular function, and therefore easy and numerically reliable to determine. By construction, if α PE = S m , we are able to predict that the solution is almost everywhere smooth and if α PE < S m that the solution possesses at least one discontinuity or a steep gradient.

Let us now illustrate the threshold α PE through an example. We represent in Fig. 7 a discontinuous solution and in Fig. 8 a smooth solution for a given mesh refinement criterion S.

• For the discontinuous solutions, as mentioned before, the threshold α = S m is not small enough to detect the region where the solution is smooth. However, if we set α = α PE < S m as a mesh refinement threshold, one can detect this region (see Fig. 7(b)).

• For smooth solutions, as emphasised in Sect. 3.3 Numerical illustration using α PE as an automatic mesh refinement threshold.

In this section we show the efficiency of the proposed method in the selection of an automatic threshold. To this end, we consider the Riemann problem (1)-( 2) with ∀x ∈ [0, 80], w(x, 0) = (h(x, 0), u(x, 0))= (5.64, 8) if x < 20 , (0.6, 8) if x > 20 . [START_REF] Kallinderis | Adaptation methods for a new navier-stokes algorithm[END_REF] for which the exact solution is composed of a left-going rarefaction and a right-going shock wave (see, for instance, Fig. 9(a), black line). Each numerical experiment starts initially with 100 cells, the maximum level of refinement is set to l max = 3, γ = 2, M = 1000 (see Eq. ( 9)) and we prescribe free boundary conditions.

The objective is to show that, for any given mesh refinement criterion, the threshold α PE is selected automatically to adapt the mesh without hand-calibration. In particular, we will confront our method to the mean method (i.e. α = β S m ) for several values of β . We recall that in this method, the threshold has to be calibrated. In this numerical experiment, we have used the following strategy: for each cell C k :

• if S n k > α, the mesh is refined and split, • if S n k < α the mesh is coarsened. which is known to be a source of spurious oscillations. We will see that our method leads to non-oscillating numerical solution without hand-calibration. For the mean method, as said before, if the threshold parameter β is not well-calibrated, then the numerical solution can yield to spurious oscillations whatever the mesh refinement criterion is. To be non redundant, we show the results only for the criterion [START_REF] Daubechies | Ten lectures on wavelets[END_REF].

In Fig. 9, we have displayed the numerical water height at time t = 2 s using the method of automatic selection of the threshold with the mesh refinement criteria (4), ( 6) and [START_REF] Ersoy | A simple and efficient new algorithm to increase the regularity and vanishing moments of biorthogonal wavelets[END_REF]. With our method, the threshold parameter α PE , being adapted at each time step, allows to detect the pertinent region of the domain to refine as one can observe in Figs. 9 and 11. It is important to note that, in spite of the mesh adaptation strategy (refined or coarsen), the numerical solution at time t = 2 s does not contain spurious oscillations. This also true at different times, as shown in Figs. 10(a)-10(d), for instance for the mesh refinement criterion (4). Furthermore, it is also interesting to highlight that, whatever the mesh refinement criterion, we obtain almost the same numerical solution with almost the same order of accuracy (as observed in Fig. 9 and Table 1). In Fig. 12, we have displayed the numerical water height at time t = 2 s using the mean method with the mesh refinement criterion [START_REF] Daubechies | Ten lectures on wavelets[END_REF]. Obviously, if the threshold parameter β is not well-chosen, we get:

Criterion N m (l) → h -h ex 1 N m (l) → u -u ex
• in Fig. 12(a), only the shock region is refined since β is set too large, • in Fig. 12(b), the numerical solution oscillates, • in Fig. 12(c), β is well-chosen, the pertinent regions to refine are well-detected and the numerical solution is well-computed (as in Fig. 9(b) with α PE ), • in Fig. 12(d), the domain is refined almost everywhere since β is set too small. Fig. 12 Numerical results for the water height at time t = 2 s using the criterion (6) with the mean method. N m is the mean number of cells used during the simulation.

Conclusion

In this paper, we have introduced a new method to set the threshold parameter automatically to get an h-AMR algorithm almost parameterless. This method allows, in principle, to refine without hand-calibration the pertinent regions. The automatic selection is constructed from a weighted function α → f (α) = αd(α) where d is decreasing rearrangement function of the mesh refinement criterion. The function f allows to compute an automatic threshold without an extra computational cost compared to the method by Dannehoffer or by Powell. Obviously other weighted functions are possible but it is important to keep in mind that a perfect automatic thresholding is a priori unreachable without any knowledge of the link between the real error and the used mesh refinement criterion. Therefore we can only focus to construct a method which has a good behaviour depending of the criterion function (with heterogeneous or homogeneous amplitudes, with smooth or discontinuous solutions, etc.).

Moreover, the method works well for any mesh refinement criterion and yields to a numerical solution without spurious oscillations. Whatever the mesh refinement criterion is, we obtain almost the same numerical solution with almost the same order of accuracy. In particular, combined with a simple gradient method as a mesh refinement criterion which is less expensive than other criteria, the h-AMR algorithm can be very efficient.

In this part, we have illustrated the efficiency on a one dimensional Riemann problem to show that our method is able to detect simple waves. In the second part of this paper [START_REF] Pons | Adaptive mesh refinement method. Part 2: Application to tsunamis propagation[END_REF], we show through several multi-dimensional test cases, including multiple waves interactions, the efficiency of the method.

Fig. 1

 1 Fig. 1 Example of hierarchical dyadic tree.

10 Fig. 2

 102 Fig. 2 Comparison of the mean value S m for two given criteria S(x) = a exp(-b(xx B ) 2 ) + exp(-5(xx A ) 2 ) where x A = 1.25 and x B = 3.75. S represents either a discontinuous solution (left) or a smooth solution (right).

Fig. 3

 3 Fig.3Illustration of the wavelet transformation for a given mesh refinement criterion computed with the Daubechies wavelets with four vanishing moments (warm colors correspond to large coefficients and cold colors to small coefficients).

Fig. 4

 4 Fig.4Illustration of the local maxima method for a mesh refinement criterion involving multiple regions.

Fig. 5

 5 Fig. 5 Dannenhoffer cumulative distribution function for a given mesh refinement criterion.

3 . 1 ,Fig. 7 Fig. 8

 3178 Fig.7The function f for the mesh refinement criterion S(x) = 200 exp(-1000(x -3.75) 2 ) + 1.25 exp(-5(x -1.25) 2 ) representing a shock type solution.

Fig. 9 Fig. 10 Fig. 11

 91011 Fig. 9 Numerical results for the water height at time t = 2 s. N m is the mean number of cells used during the simulation.

  x) on adaptive mesh with l max = 3 S k n /S m with β = 0.00000001 and N m = 270 Refined zones (d) β = 10 -7

Table 1

 1 Convergence rate of the L 1 discretisation error obtained with the refinement level l = 1, 2, 3 and 4.

	1

the wavelet coefficients in this example are computed with the Daubechies wavelet with four vanishing moments[START_REF] Daubechies | Ten lectures on wavelets[END_REF][START_REF] Ersoy | A simple and efficient new algorithm to increase the regularity and vanishing moments of biorthogonal wavelets[END_REF]

Acknowledgements This work is partially supported the Project MTM2011-29306-C01-01 from the MICINN (Spain) and the French national research project TANDEM (Tsunamis in the Atlantic and the English channel : definition of the effects through numerical modeling), the French government (Projets Investissement d'Avenir, agreement reference number ANR-11-RSNR-0023-01). The authors wish to thank J-J. Alibert for its valuable remarks and comments. The author also wish to thanks the referees for their valuable remarks which led to substantial improvement of the first version of this paper.


Afzal Khan

Viet Huong Nguyen

David Muñoz-Rojas
email: david.munoz-rojas@grenoble-inp.fr

Sara Aghazadehchors

Carmen Jiménez

Ngoc Duy

Daniel Bellet

Stability enhancement of silver nanowire networks with conformal ZnO coatings deposited by atmospheric pressure spatial atomic layer deposition

Keywords: Transparent electrode, stability, conformal coating, ZnO, Spatial Atomic Layer Deposition, AP-SALD, failure, silver nanowire

Silver nanowire (AgNW) networks offer excellent electrical and optical properties and have emerged as one of the most attractive alternatives to transparent conductive oxides to be used in flexible optoelectronic applications. However, AgNW networks still suffer from chemical, thermal and electrical instabilities which in some cases can hinder their efficient integration as transparent electrodes in devices such as solar cells, transparent heaters, touch screens or organic light emitting diodes (OLEDs). We have used atmospheric pressure spatial atomic layer deposition (AP-SALD) to fabricate hybrid transparent electrode materials in which the AgNW network is protected by a conformal thin zinc oxide layer. The choice of AP-SALD allows to maintain the low-cost and scalable processing of AgNW based transparent electrodes. The effects of the ZnO coating thickness on the physical properties of AgNW networks are presented. The composite electrodes show a drastic enhancement of both thermal and electrical stabilities. We found that bare AgNWs were stable only up to 300 °C when subjected to thermal ramps while the ZnO coating improved stability up to 500 °C.

Similarly, ZnO coated AgNWs exhibited an increase of a 100 % in electrical stability with respect to bare networks, withstanding up to 18 V. A simple physical model shows that the origin of the stability improvement is the result of hindered silver atomic diffusion thanks to the presence of the thin oxide layer and the quality of the interfaces of hybrid electrodes. The effects of ZnO coating on both the network adhesion and optical transparency are also discussed. Finally, we show that the AP-SALD ZnO-coated AgNW networks can be effectively used as very stable transparent heaters.

Introduction

Transparent electrodes (TE) constitute a critical component within a wide range of devices including smart windows, touch screens, solar cells, organic light emitting diodes (OLED) or transparent heaters. Such electrodes need to be optically transparent and electrically conductive to a certain extent [START_REF] Ellmer | Past Achievements and Future Challenges in the Development of Optically Transparent Electrodes[END_REF][START_REF] Granqvist | Transparent Conductors as Solar Energy Materials: A Panoramic Review[END_REF] depending on the targeted application. While the main transparent conductive materials studied for several decades have been transparent conductive oxides (TCO), alternative emerging thin layers have been explored recently. [START_REF] Hecht | Emerging Transparent Electrodes Based on Thin Films of Carbon Nanotubes, Graphene, and Metallic Nanostructures[END_REF] Among them, metallic nanowire networks and specially those composed of silver nanowires (AgNW) have attracted much attention in the past few years. [START_REF] De | Silver Nanowire Networks as Flexible, Transparent, Conducting Films: Extremely High DC to Optical Conductivity Ratios[END_REF][START_REF] Langley | Flexible Transparent Conductive Materials Based on Silver Nanowire Networks: A Review[END_REF][START_REF] Guo | Flexible Transparent Conductors Based on Metal Nanowire Networks[END_REF][START_REF] Sannicolo | Metallic Nanowire-Based Transparent Electrodes for Next Generation Flexible Devices : A Review[END_REF][START_REF] He | Nanostructured Transparent Conductive Films: Fabrication, Characterization and Applications[END_REF] By virtue of their excellent electrical and optical properties (10-20 Ω/sq sheet resistance and about 90 % transmittance in the visible range) [START_REF] De | Silver Nanowire Networks as Flexible, Transparent, Conducting Films: Extremely High DC to Optical Conductivity Ratios[END_REF][START_REF] Lagrange | Optimization of Silver Nanowire-Based Transparent Electrodes: Effects of Density, Size and Thermal Annealing[END_REF] and their mechanical flexibility, [START_REF] Sannicolo | Metallic Nanowire-Based Transparent Electrodes for Next Generation Flexible Devices : A Review[END_REF] silver nanowire networks are considered as promising candidates for flexible optoelectronic applications. Another clear advantage of AgNW based TE is the possibility to establish simple large-scale fabrication through solution processing [START_REF] Bellet | Transparent Electrodes Based on Silver Nanowire Networks: From Physical Considerations towards Device Integration[END_REF] and a lower price than indium tin oxide (ITO), [START_REF] Emmott | Environmental and Economic Assessment of ITO-Free Electrodes for Organic Solar Cells[END_REF] which is so far considered as the most standard and efficient transparent conducting oxide. [START_REF] Ellmer | Past Achievements and Future Challenges in the Development of Optically Transparent Electrodes[END_REF] It has even been shown that by further optimizing the deposition and post-deposition conditions, and the aspect ratio (length/diameter) of the nanowires, the resulting networks can surpass the electro-optical properties of (ITO). [START_REF] Lagrange | Optimization of Silver Nanowire-Based Transparent Electrodes: Effects of Density, Size and Thermal Annealing[END_REF] AgNW networks have been investigated or even efficiently integrated into devices such as transparent heaters, [START_REF] Celle | Highly Flexible Transparent Film Heaters Based on Random Networks of Silver Nanowires[END_REF][START_REF] Kiruthika | Large Area Defrosting Windows Based on Electrothermal Heating of Highly Conducting and Transmitting Ag Wire Mesh[END_REF][START_REF] Sorel | Relationship between Material Properties and Transparent Heater Performance for Both Bulk-like and Percolative Nanostructured Networks[END_REF][START_REF] Ergun | High-Performance, Bare Silver Nanowire Network Transparent Heaters[END_REF] solar cells, [START_REF] Guo | ITO-Free and Fully Solution-Processed Semitransparent Organic Solar Cells with High Fill Factors[END_REF][START_REF] Guo | Fully Printed Organic Tandem Solar Cells Using Solution-Processed Silver Nanowires and Opaque Silver as Charge Collecting Electrodes[END_REF][START_REF] Yu | Silver Nanowire-Polymer Composite Electrodes for Efficient Polymer Solar Cells[END_REF][START_REF] Morgenstern | Ag-Nanowire Films Coated with ZnO Nanoparticles as a Transparent Electrode for Solar Cells[END_REF][START_REF] Langley | Silver Nanowire Networks: Physical Properties and Potential Integration in Solar Cells[END_REF] touch screens, [START_REF] Lee | Very Long Ag Nanowire Synthesis and Its Application in a Highly Transparent, Conductive and Flexible Metal Electrode Touch Panel[END_REF][START_REF] Wang | Silver Nanowire Electrodes: Conductivity Improvement Without Post-Treatment and Application in Capacitive Pressure Sensors[END_REF] electromagnetic shielding [START_REF] Kim | Transparent and Flexible Film for Shielding Electromagnetic Interference[END_REF] or antennas, [START_REF] Song | Stretchable and Reversibly Deformable Radio Frequency Antennas Based on Silver Nanowires[END_REF] OLEDs, [START_REF] Lee | Highly Efficient and Low Voltage Silver Nanowire-Based OLEDs Employing a n-Type Hole Injection Layer[END_REF] in electrochromic devices, [START_REF] Yuksel | A New High-Performance Blue to Transmissive Electrochromic Material and Use of Silver Nanowire Network Electrodes as Substrates[END_REF] for medical purposes [START_REF] Choi | Stretchable Heater Using Ligand-Exchanged Silver Nanowire Nanocomposite for Wearable Articular Thermotherapy[END_REF] as well electrically conductive textiles [START_REF] Atwa | Silver Nanowire Coated Threads for Electrically Conductive Textiles[END_REF] . Despite the attracting potentialities described here above, the integration of AgNW based TE in real devices is not yet widely considered in industrial devices due namely to their potential thermal [START_REF] Langley | Metallic Nanowire Networks: Effects of Thermal Annealing on Electrical Resistance[END_REF] and electrical [START_REF] Khaligh | Failure of Silver Nanowire Transparent Electrodes under Current Flow[END_REF] instabilities, as well as low adhesion [START_REF] Sannicolo | Metallic Nanowire-Based Transparent Electrodes for Next Generation Flexible Devices : A Review[END_REF] and aging issues. [START_REF] Mayousse | Stability of Silver Nanowire Based Electrodes under Environmental and Electrical Stresses[END_REF][START_REF] Celle | Oxidation of Copper Nanowires Based Transparent Electrodes in Ambient Conditions and Their Stabilization by Encapsulation. Application to Transparent Film Heaters[END_REF] For instance, while an electrical welding [START_REF]Joining of Individual Silver Nanowires via Electrical Current[END_REF] or thermal annealing [START_REF] Langley | Metallic Nanowire Networks: Effects of Thermal Annealing on Electrical Resistance[END_REF][START_REF] Sannicolo | Direct Imaging of the Onset of Electrical Conduction in Silver Nanowire Networks by Infrared Thermography: Evidence of Geometrical Quantized Percolation[END_REF] is beneficial for reducing junction resistance (between adjacent AgNWs) thanks to local sintering, it is also known that reaching too high a temperature leads to spheroidization of the AgNWs. This causes the loss of the percolating nature of the network, thus eventually leading to an infinite electrical resistance. [START_REF] Langley | Metallic Nanowire Networks: Effects of Thermal Annealing on Electrical Resistance[END_REF] Left unprotected, AgNW networks would undergo local oxidation and/or melt, which adversely affect the conductivity of AgNW films. While the melting point of bulk silver is above 960 °C, degradation induced by thermal annealing can occur at much lower temperature (i.e. < 300 °C) due to the high surface-to-volume ratio in nanowires. [START_REF] Lagrange | Optimization of Silver Nanowire-Based Transparent Electrodes: Effects of Density, Size and Thermal Annealing[END_REF] The idea of preparing composite electrodes by coating metallic nanowires with a very thin oxide layer has already been considered in literature [START_REF] Morgenstern | Ag-Nanowire Films Coated with ZnO Nanoparticles as a Transparent Electrode for Solar Cells[END_REF] as a possible remedy to such instabilities. Some researchers have used titanium dioxide (TiO 2 ) coatings [START_REF] Ramasamy | Effects of TiO2 Shells on Optical and Thermal Properties of Silver Nanowires[END_REF][START_REF] Song | Highly Robust Silver Nanowire Network for Transparent Electrode[END_REF][START_REF] Zhu | Fused Silver Nanowires with Metal Oxide Nanoparticles and Organic Polymers for Highly Transparent Conductors[END_REF] , as well as ZnO [START_REF] Kim | Highly Transparent Low Resistance ZnO / Ag Nanowire / ZnO Composite Electrode for Thin Film Solar Cells[END_REF][START_REF] Chen | Thermally Stable Silver Nanowire-Polyimide Transparent Electrode Based on Atomic Layer Deposition of Zinc Oxide on Silver Nanowires[END_REF][START_REF] Pham | Enhancement of the Electrical Properties of Silver Nanowire Transparent Conductive Electrodes by Atomic Layer Deposition Coating with Zinc Oxide[END_REF] or aluminum doped zinc oxide (ZnO:Al). [START_REF] Göbelt | Encapsulation of Silver Nanowire Networks by Atomic Layer Deposition for Indium-Free Transparent Electrodes[END_REF] Clearly, the shell used has shown encouraging results in terms of stability enhancement. However, in order for such composite electrodes to remain attractive in terms of cost and simple processing, the oxide layer should be deposited by low-cost, vacuum-free scalable approaches. The methods so far reported (Atomic Layer Deposition (ALD), spin-coating or sputtering) [START_REF] Ramasamy | Effects of TiO2 Shells on Optical and Thermal Properties of Silver Nanowires[END_REF][START_REF] Song | Highly Robust Silver Nanowire Network for Transparent Electrode[END_REF][START_REF] Zhu | Fused Silver Nanowires with Metal Oxide Nanoparticles and Organic Polymers for Highly Transparent Conductors[END_REF][START_REF] Kim | Highly Transparent Low Resistance ZnO / Ag Nanowire / ZnO Composite Electrode for Thin Film Solar Cells[END_REF][START_REF] Chen | Thermally Stable Silver Nanowire-Polyimide Transparent Electrode Based on Atomic Layer Deposition of Zinc Oxide on Silver Nanowires[END_REF][START_REF] Göbelt | Encapsulation of Silver Nanowire Networks by Atomic Layer Deposition for Indium-Free Transparent Electrodes[END_REF] can be associated to drawbacks such as the use of vacuum or the difficulty for upscaling, thus compromising the intrinsic low-cost and scalability of AgNW based TE. Atmospheric Pressure Spatial Atomic Layer Deposition (AP-SALD) is a promising recent approach to ALD that can operate at atmospheric pressure, is easily scalable and up to two orders of magnitude faster than ALD. [START_REF] Muñoz-Rojas | Spatial Atmospheric Atomic Layer Deposition: A New Laboratory and Industrial Tool for Low-Cost Photovoltaics[END_REF] This is so thanks to the spatial separation of precursors in AP-SALD, which eliminates the need for purges between successive precursor injections in conventional ALD. [START_REF] Muñoz-Rojas | Spatial Atmospheric Atomic Layer Deposition: A New Laboratory and Industrial Tool for Low-Cost Photovoltaics[END_REF][START_REF] Illiberi | Spatial Atomic Layer Deposition of Zinc Oxide Thin Films[END_REF][START_REF] Poodt | Spatial Atomic Layer Deposition: A Route towards Further Industrialization of Atomic Layer Deposition[END_REF] AP-SALD has already been used to deposit high-quality oxide films as passive and active components for new generation solar cells and other devices. [START_REF] Muñoz-Rojas | Spatial Atmospheric Atomic Layer Deposition: A New Laboratory and Industrial Tool for Low-Cost Photovoltaics[END_REF][START_REF] Hoye | Research Update: Atmospheric Pressure Spatial Atomic Layer Deposition of ZnO Thin Films: Reactors, Doping, and Devices[END_REF] The efficiency of AP-SALD films to protect transparent heaters based on Cu nanowires from chemical oxidation has also been shown recently. [START_REF] Celle | Oxidation of Copper Nanowire Based Transparent Electrodes in Ambient Conditions and Their Stabilization by Encapsulation : Application to Transparent Fi Lm Heaters[END_REF] The present study focuses on the use of AP-SALD to deposit ZnO coatings on AgNWs based TE, and the thorough investigation of the effects of the ZnO coating thickness on the properties of AgNW networks in terms of electrical and thermal stabilities. We have used a home-made AP-SALD system to fabricate ZnO-AgNW composite electrodes. We demonstrate that AP-SALD is a suitable technique to deposit in a few minutes a conformal compact oxide layer of a controlled thickness (up to 30 nm), which drastically improves the thermal and electrical stability of the networks. We also introduce a simple model that allows to explain the stability enhancement provided by the ZnO coating. The effects of ZnO coating on both the network adhesion and optical transparency are also discussed. Finally, we show that the AP-SALD ZnO-coated AgNW networks can be effectively used as very stable transparent heaters.

Methods and experimental conditions

Deposition of AgNW networks:

AgNW suspensions in isopropanol were purchased from Seashell technology, USA. The average length of AgNWs was 32 µm while their average diameter was 130 nm. The asreceived suspensions, with initial concentration of 10 mg/ml, were diluted to various concentrations and then deposited on Corning glass (C1737-S111) substrates using spin coating at room temperature. The substrates were cleaned through a sonication in acetone, then rinsed in isopropanol, then washed with distilled water and finally dried with N2 gas before AgNW deposition. Deposition parameters (nanowires solution concentration, rotating speed of the spin coater and dispensing rate) and post-deposition thermal annealing parameters (temperature and duration) were optimized to obtain highly transparent (~90 % transmittance) and highly conductive (~10 Ω/sq sheet resistance) AgNW networks. The best electro-optical properties of AgNW networks were obtained by spin-coating 1 ml of suspension (0.5 mg/ml, in two stages of 0.5 ml each with a 20 s pause between them, allowing the surface to dry) over 2.5 cm x 2.5 cm corning glass substrates at ambient atmospheric conditions using a spin-coater with rotation of 1500 rpm. The optimal postdeposition thermal annealing (250 °C for 30 minutes) step is necessary to desorb the organic residues and to decrease the electrical resistance of the junctions between adjacent AgNWs. [START_REF] Langley | Metallic Nanowire Networks: Effects of Thermal Annealing on Electrical Resistance[END_REF] A second set of samples was deposited for performing isothermal annealing tests and to test the electrical stability using different voltage ramps. In this case, the diameter of NWs was around 90 nm and the lengths could vary from 10 to 20 um. The initial suspension was diluted to 0.2g/L prior to deposition. An airbrush infinity CR was used to spay NWs on corning glass substrates. During deposition, a N 2 flow at 3 bar was used, and the airbrush scanned the substrate in the X-Y direction for 15 cycles. The substrates were kept at 90 °C during the deposition. The as-deposited samples were annealed at 230 °C for 1 h to minimize the electrical resistance and ensure intimate junction between adjacent AgNWs.

Deposition of ZnO coatings

The ZnO coatings were deposited using a home-made AP-SALD system previously reported. [START_REF] Muñoz-Rojas | Spatial Atmospheric Atomic Layer Deposition: A New Laboratory and Industrial Tool for Low-Cost Photovoltaics[END_REF][START_REF] Nguyen | Deposition of ZnO Based Thin Films by Atmospheric Pressure Spatial Atomic Layer Deposition for Application in Solar Cells[END_REF][START_REF] Muñoz-Rojas | Spatial Atomic Layer Deposition (SALD), an Emerging Tool for Energy Materials[END_REF] The deposition parameters were optimized and then ZnO growth was carried out at 200 °C, using diethylzinc ((C 2 H 5 ) 2 Zn; DEZ) and water vapor as precursors for zinc and oxygen, respectively. The samples were placed at a distance of 200 µm under the injection head oscillating at 10 cm/s. In our deposition head, there are two channels for the metal precursor. Thus, each sample oscillation (forth and back) corresponds to four ALD cycles. [START_REF] Muñoz-Rojas | High-Speed Atmospheric Atomic Layer Deposition of Ultra Thin Amorphous TiO 2 Blocking Layers at 100 °C for Inverted Bulk Heterojunction Solar Cells[END_REF] The system allows to deposit over an area of 5 cm x 5 cm.

Morphological and optical characterization

Optical characterization of the AgNW networks was performed with a Perkin Elmer Lambda 950 UV-Visible-Near IR spectrophotometer equipped with an integrating sphere, which allows measuring the total optical transmittance/reflection. Such measurements were performed in two different locations for each sample and then averaged.

Optical microscopy images of the AgNW networks were obtained with a Leica DM LM microscope. Scanning electron microscopy (SEM) and transmission electron microscopy (TEM) images were obtained using, respectively, a FEG-SEM Environmental FEI QUANTA 250 operating at an accelerating voltage of 10 keV and a JEOL JEM-2010 microscope operating at 200 keV. The sample preparation for TEM imaging was performed using the crushing technique.

Thermal and electrical stability studies

The stability of AgNW networks was assessed through in situ measurement of the network electrical resistance during either thermal or voltage ramps by using a two-probe setup with a Keithley 2400 sourcemeter. Thermal stability was studied by subjecting the AgNW networks to thermal ramps from room temperature to 500 o C at a rate of 2 o C/min. Isothermal test at 320 °C for 1 hour w also performed.

The electrical stability of the networks was measured in two ways. In the first constant voltage ramp at a rate of 0.1 volt/min, were applied up to specimen failure (i.e. when the electrical resistance diverged) to samples coated with different ZnO thicknesses (the networks used were fabricated using the AgNWs with an average diameter of 130 nm). In the second, equivalent networks made using the AgNWs with average diameter of 90 nm where coated with 25 nm of ZnO and different voltage ramps were applied, namely, 0.05, 0.075 and 0.1 volt/min.

To assess the use of these AgNW networks as transparent heater, different voltages were applied while the temperature of the samples was measured in situ using a thin resistance temperature detector (RTD) embedded in a 120 µm thick PET layer, therefore limiting thermal inertia. The RTD was placed in thermal contact with the sample, below the substrate (details on the setup have been reported previously. [START_REF] Lagrange | Optimization of Silver Nanowire-Based Transparent Electrodes: Effects of Density, Size and Thermal Annealing[END_REF][START_REF] Lagrange | Understanding the Mechanisms Leading to Failure in Metallic Nanowire-Based Transparent Heaters, and Solution for Stability Enhancement[END_REF]

Results and Discussion

ZnO coating of AgNW network with AP-SALD and effect on adhesion and transparence

A SEM image of the bare AgNWs annealed at 250 o C for 30 minutes is shown in Figure 1a, where effective sintering of the wires at the junctions is clearly visible. A coated AgNW network with a 30 nm thick ZnO coating deposited with AP-SALD is shown in the SEM image in Figure 1b. The distinct morphology of the ZnO coating can be observed both on the nanowires and the glass surface. As shown in the figures, the surface of the bare nanowires appears smooth while that of ZnO-coated AgNWs show the typical granular morphology of ZnO coatings deposited by AP-SALD. A transmission electron microscopy (TEM) image of a 30 nm coating of ZnO on a AgNW is shown in Figure 1c. Such coating appears very conformal: this point is crucial since this is a prerequisite for an efficient protection against any instability. As in ALD, in AP-SALD the ZnO layer thickness is proportional to the number of AP-SALD cycles, i.e. sample oscillations under the head. The growth rate r was thus estimated from the thickness measured for the thicker ZnO layer deposited on AgNWs divided by the number of AP-SALD oscillations. As reported in Figure 1c, a ZnO thickness of 30 nm is observed for 60 AP-SALD cycles. The growth rate was thus ~ 0.50 nm/oscillation (as stated above, in our AP-SALD system a sample oscillation is equivalent to 4 conventional ALD cycles [START_REF] Nguyen | Deposition of ZnO Based Thin Films by Atmospheric Pressure Spatial Atomic Layer Deposition for Application in Solar Cells[END_REF] , thus yielding a growth per cycle (GPC) of about 1.25 Å/ALD cycle). Figure 1d shows a back scattered SEM image of a ZnO coated AgNW which enables to highlight the elemental contrast: the usual fivefold symmetry twinning of AgNWs 7 is well seen, as well as the surrounding uniform ZnO layer. Again the conformal and continuous aspect of the coating appears in a very clear way in this image.

Before focusing on the effects of ZnO coating on electrical and thermal stability of AgNW networks, we first consider the effects on both the network adhesion and optical transparency.

It is known that low adhesion of AgNW networks to the substrate can be one of the drawbacks of these materials when used as TE [START_REF] Sannicolo | Metallic Nanowire-Based Transparent Electrodes for Next Generation Flexible Devices : A Review[END_REF] : even a gentle rub can remove nanowires from the substrate. The influence of the AP-SALD ZnO coating on the adherence of AgNW network on the glass substrate was thus investigated. To do so a scotch tape test was applied to both the bare and the coated networks to qualitatively compare their adhesion strength to the substrates. Figure 2a and2b show the optical microscope images of the bare and the ZnO coated AgNWs networks at the edge of the locations where the tape was applied. The upper parts in the figures correspond to the areas in which the scotch tape test was applied, while the lower parts correspond to areas outside the tested zone. Figure 2a The effect of the ZnO coating on electrode transparency was also evaluated. Figure 2c reports the influence of the ZnO coating thickness on the optical transparency (measured for a wavelength of 550 nm) of coated AgNW networks. Note the network transmittance is measured after subtracting the substrate contribution, as explained in Lagrange's et al.

article. [START_REF] Lagrange | Optimization of Silver Nanowire-Based Transparent Electrodes: Effects of Density, Size and Thermal Annealing[END_REF] In a first approximation a linear decrease of optical transparency from 93.5 % to 82.5 % can be observed (as expected for a linear expansion of an exponential function close to unity). For an optimized integration of coated AgNW network a compromise should be considered, depending on the application targeted. Indeed, as seen below, a thicker ZnO coating is associated to a more efficient thermal and electrical stability but at the expenses of optical transparency. This slight decrease of optical transparency is nevertheless still acceptable for many applications. Finally it is worth noticing that the electrical resistance of AgNW network remains similar, within experimental incertitude, before and after ZnO coating.

Effect of the ZnO coating on the electrical and thermal stability of AgNW networks

For all considered applications, both electrical and thermal instabilities of AgNW networks can be a key issue for their efficient integration into real devices, since this could drastically reduce the lifetime of the TE and thus of the associated device. The effect of ZnO coating on electrical and thermal stabilities of the AgNW networks has been evaluated as a function of the thickness of the ZnO coating. In this regard, there is no universal definition of metallic network stability. In the present paper we define the stability of an AgNW network (during either thermal or voltage ramps) as a domain of variation of voltage or temperature where the change of network resistance is reversible. Inversely the instability can be observed when the electrical resistance suddenly increases during a ramp from few tens of Ohms to more than 10 4 Ohms. We first consider thermal stability and then the electrical stability.

Similar AgNW networks were spin-coated on Corning glass substrates, then ZnO coatings of various thicknesses were deposited by AP-SALD. The initial electrical resistance of each specimen was similar with an average value of 12.1  and a standard deviation of 2.1  regardless the ZnO layer thickness. The stability of AgNW networks to high temperatures was investigated during thermal ramps from room temperature to 500 o C at a heating rate of 2 o C/min. Figure 3a shows that in all cases AgNW networks exhibit a slight linear increase of resistance with temperature, which is related to the electron-phonon interaction [START_REF] Bid | Temperature Dependence of the Resistance of Metallic Nanowires of Diameter ⩾15nm: Applicability of Bloch-Grüneisen Theorem[END_REF] . For bare AgNW networks and those coated with a thin ZnO layer, this linear dependence is then followed by a very sharp increase of the electrical resistance. This last feature stems from a non-reversible phenomenon which leads to the spheroidization of the AgNWs (due to the total surface energy reduction). This morphological instability is called the Plateau-Rayleigh instability and has been already well described in the case of metallic nanowires [START_REF] Lagrange | Optimization of Silver Nanowire-Based Transparent Electrodes: Effects of Density, Size and Thermal Annealing[END_REF][START_REF] Langley | Metallic Nanowire Networks: Effects of Thermal Annealing on Electrical Resistance[END_REF][START_REF] Urban | Interplay of Rayleigh and Peierls Instabilities in Metallic Nanowires[END_REF][START_REF] Karim | Morphological Evolution of Au Nanowires Controlled by Rayleigh Instability[END_REF] . The spheroidization is at the origin of the loss of the electrical percolating nature of the network, which can then be observed through the sharp divergence in electrical resistance for some specimens. Figure 3a clearly shows that by increasing the thickness the ZnO coating, higher temperatures can be achieved before the instability is observed. For instance, a bare AgNW network deteriorates at 315 °C while the AgNW network coated with ZnO of at least 25 nm can sustain temperatures up to 500 o C while showing finite electrical resistance. In particular, networks coated with only 15 nm show a resistance that is 5 times the original one after undergoing the thermal ramp. With 20 nm of ZnO coating, the increase of resistance goes down to 1.8 times the original value. Finally, for 25 nm and thicker coatings, the final resistance after the thermal ramp is only 1.1 times higher than the initial value. The effect of the thickness of the ZnO coating appears quite evident in terms of resistance reversibility. In order to assess the network thermal stability by another method, the resistance of the networks was measured versus time at constant temperature (320 °C). This feature can be observed in Figure 3c, where the relative variation of electrical resistance R/R 0 is plotted against time during an annealing at 320 °C for specimens associated to several ZnO coating thicknesses. Clearly, increasing the ZnO coating enhances the stability of the network. A ZnO layer thickness of 20-30 nm appears sufficient to maintain a constant electrical resistance.

In a similar approach, we subjected the networks to electrical stress to study their electrical stability. Figure 4a reports the variation of the electrical resistance of bare and coated (different thicknesses) AgNW network versus applied voltage during a voltage ramp of 0.1 V/min. In all cases a slight increase of the electrical resistance versus voltage is observed at low voltage. The origin of such a resistance increase is the Joule effect: when an electric current flows through the network the temperature of the latter rises and the electrical resistance as well, due to the electron-phonon interactions. This evolution is reversible, at least when low voltages are considered. For larger voltages, the network can exhibit an electrical failure which is associated with a sudden and non-reversible increase of the electrical resistance. The resistance curves obtained during the in situ electrical measurements displayed in Figure 4a shows that again increasing the thickness of the ZnO coating allows to reach higher voltages before failure takes place... This is similar to what has been observed above when considering the thermal stability. While bare networks show failure at around 9 volts, the stability of ZnO coated networks reaches up to 18 V for 30 nm of ZnO coating, which corresponds to a 100% increase of the associated voltage failure. Figure 4b In one of our previous studies [START_REF] Lagrange | Understanding the Mechanisms Leading to Failure in Metallic Nanowire-Based Transparent Heaters, and Solution for Stability Enhancement[END_REF] we reported a physical model estimating the link between the average electrical current density and the voltage V applied between two opposite electrodes separated by a distance W. The relation can be written as follows:

where  is the electrical resistivity of the network. Therefore from the voltage values observed in Figure 4 one can deduce that the current density at failure varies between 1.6 10 10

A/m 2 for non-coated AgNW and 3. AgNW was clearly observed due to electromigration. [START_REF] Stahlmecke | Electromigration in Self-Organized Single-Crystalline Silver Nanowires[END_REF] The fact that the current density associated with the network failure is lower than the one observed by Stahlmecke et al. might indicate that the AgNW associated with smaller diameters degrade first, or areas of the networks with the lower-than-average nanowire density.

Rationalization of the effect of the ZnO coating on the stability of AgNW networks

The origin of the nanowire Plateau-Rayleigh instability has already been reported several times, both for metallic nanowires and as well for other types of nanowires. [START_REF] Urban | Interplay of Rayleigh and Peierls Instabilities in Metallic Nanowires[END_REF][START_REF] Karim | Morphological Evolution of Au Nanowires Controlled by Rayleigh Instability[END_REF][START_REF] Xue | Engineering Island-Chain Silicon Nanowires via a Droplet Mediated Plateau-Rayleigh Transformation[END_REF] The main driving force at the origin of such instability is the reduction in total surface energy. Figure 5a shows a simple schematic representation of a junction between two AgNWs. When subjected to a thermal annealing or a voltage ramp (the Joule effect increasing the temperature of the network), the junction undergoes a morphological evolution which aims at decreasing its total surface energy. The consequence is local sintering, and therefore the junction resistance is lowered, as well as the network electrical resistance [START_REF] Langley | Metallic Nanowire Networks: Effects of Thermal Annealing on Electrical Resistance[END_REF] . Figure 5b and 5c illustrate such morphological changes at different stages of the treatment: the first one leads to the minimum electrical resistance (local sintering) of the junction while the second one is associated to the onset of network degradation (with local junction breakage probably at the hottest points of the AgNW network). An example of actual morphological instability is presented in Figure 5d, which shows a SEM picture of a AgNW network subjected to an electrical load or to a thermal annealing close to the failure limit (before a complete spheroidization). All these morphological modifications are possible provided atomic surface diffusion can take place.

But such atomic diffusion can be hindered by an oxide coating. This is the reason, for instance, why sintering of metallic powder is operated in a reductive atmosphere (often with presence of hydrogen), [START_REF] Vagnon | 3D Statistical Analysis of a Copper Powder Sintering Observed in Situ by Synchrotron Microtomography[END_REF][START_REF] Vagnon | Deformation of Steel Powder Compacts during Sintering: Correlation between Macroscopic Measurement and in Situ Microtomography Analysis[END_REF] to prevent the formation of surface oxide layer that would prevent the sintering. A thin interfacial SiO 2 layer is also known to weaken the thermal diffusion of metallic boron or arsenic into silicon. [START_REF] Nguyen | Use of P-and n-Type Vapor Phase Doping and Sub-Melt Laser Anneal for Extension Junctions in Sub-32 Nm CMOS Technology[END_REF][START_REF] Takeuchi | Vapor Phase Doping with N-Type Dopant into Silicon by Atmospheric Pressure Chemical Vapor Deposition[END_REF] where is the voltage increase ramp. Since the ZnO coating induces an increase of the failure time (see Figure 4b), then dependence with ZnO coating thickness (L ZnO ) can be written as follows:

eq.( 1)

where and are the failure time (for which the AgNW network electrical resistance is diverging) for bare AgNW and ZnO coated AgNW networks, respectively.

As shown schematically in Figure 5, can be associated to the time (or voltage) required for Ag atoms to undergo enough displacement enabling the morphological instability of bare AgNW. When an AgNW is coated with a thin oxide layer, the additional diffusion of Ag atoms through the thin oxide layer induces a delay, noted . By using the basic diffusion laws in first approximation, the latter is directly dependent on the oxide layer thickness ( ) through:

eq.( 2
)
where D is the diffusion coefficient of Ag atoms through the oxide barrier. Therefore this simple model would lead to a voltage failure versus ZnO coating thickness L ZnO given by: eq.( 3)

where is the voltage rate associated to the voltage ramp (equal to 0.1 V/min).

Figure 5g reports the same experimental data as already presented in Figure 4b along with a fit using eq.( 3). In spite of the simplicity of the model, the agreement appears very good between the model and experimental data. The only unknown parameter in eq.( 3) is the diffusion coefficient D of silver atoms though ZnO oxide. The fit presented in Figure 5g yields a D value of 9.9 10 -16 cm 2 /s. The comparison between this D value and the ones reported in literature is not trivial since the diffusion coefficient depends drastically on: i/ the structural properties of the oxide (presence of defects such as grain boundaries or voids) which themselves depend on the growth experimental conditions, and ii/ temperature since diffusion is a thermally activated process.

The diffusion coefficient of Ag atoms in ZnO, D, has a very low value since Ag ion has a low solubility limit in ZnO at equilibrium. This stems from the different valence states of Ag ion and Zn as well as from the larger Ag ionic radius when compared with Zn 61 . Sakaguchi et al. investigated the diffusion of Ag in polycrystalline ZnO in the temperature range 700-900 °C. [START_REF] Sakaguchi | Ion Implantation and Diffusion Behavior of Silver in Zinc Oxide[END_REF] They found a Ag diffusion coefficient in ZnO equal to 7.8 10 -13 cm 2 /s at 900 °C. For a lower temperature range, McBrayer et al. [START_REF] Mcbrayer | Diffusion of Metals in Silicon Dioxide[END_REF] studied metal diffusion in SiO2 oxide layers and reported that Ag diffusion coefficient D is equal to 10 -15 cm 2 /s at a temperature of 221 °C, which is typically the measured temperature associated to the AgNW network failure. These D values from literature are thus in agreement with the value obtained from our model (3) (see Figure 5g).

According to our model, the failure voltage should also be affected by the value of , i.e.

the rate at which the voltage is increased. Thus, in order to further test the validity of the model, we also performed voltage ramps using different values on networks coated with 25 nm of ZnO. As shown in Figure 5h), a steep increase in resistance at about 11 V is observed when using a of 0.05 V/min. For faster voltage ramps, the failure voltage increased from 12

V for a of 0.075 V/min to over 14 V for of 0.1 V/min. Figure 5i) shows the failure voltage versus the product •(L ZnO ) 2 , with data both from constant voltage ramps with different ZnO thicknesses (Figure 4a)) and from voltage ramps with different rates (Figure 5h)). Again, the agreement appears very good between the model and experimental data, even more taking into account that the networks used in both cases (Figure 4a) and Figure 5h))

were fabricated using AgNWs with different diameter and lengths, and with different deposition methods (spin-coating vs. Spray-coating), respectively, as detailed in the experimental section

Evaluation of the performance of transparent heaters based on ZnO coated

AgNWs

One of the key issues associated with AgNW networks when used for transparent heater applications is their electrical non-reversibility after applying voltage ramps or voltage plateaus. [START_REF] Lagrange | Understanding the Mechanisms Leading to Failure in Metallic Nanowire-Based Transparent Heaters, and Solution for Stability Enhancement[END_REF] To assess the reversibility of the electrical resistance of bare and ZnO coated AgNW networks, voltage plateaus were applied at 1, 3, 5, 7, 5, 3 and 1 V with a duration of 45 minutes each. Figure 6a shows the in situ measured electrical resistance of the studied AgNW networks. The bare AgNW network clearly shows a continuous degradation especially for the largest voltage plateau (7 V) with a continuous observed increase of the resistance with time. Conversely, a ZnO coated network (25 nm thick) became reversible and stable under the same conditions. Another way of demonstrating these effects is to plot the measured electrical resistance versus temperature. As shown in Figure 6b, the ZnO coating has clear positive effects on the reversibility of the electrical resistance when compared with the bare network. Interestingly, Figure 6b shows that AgNW networks, if well protected against morphological instability, can not only be used as transparent heaters but as well as a transparent temperature sensors. Thanks to the positive effects of the ZnO coating, which increases the voltage failure limits, the heating capability of AgNW network can be enhanced. One such example is shown in Figure 7, where for a voltage ramp of 0.1 V/min the bare AgNW network degraded at around 8 V and the maximum associated temperature measured was 150 o C. Conversely, the network protected by a 25 nm thick ZnO coating degraded at around 14 V and the related measured maximum temperature was 350 o C (an increase of 133%). The maximum temperature limit can therefore be drastically enhanced thanks to the ZnO coating over the AgNWs since atomic surface diffusion can be partially delayed.

Conclusions

In summary, ZnO thin layers of different thicknesses were deposited on silver nanowire (AgNW) networks by atmospheric pressure spatial atomic layer deposition (AP-SALD) and the effects on the main physical properties of the networks were studied. It is shown that the ZnO coatings drastically improve adhesion as well as thermal and electrical stability of AgNW networks, with only a small decrease in optical transparency. It is worth mentioning that, as opposed to previous works, ZnO coating has been deposited by AP-SALD, which is a new approach to ALD that maintains the high quality of coatings associated to ALD but that can be performed at atmospheric pressure, offering thus deposition rates orders of magnitude faster than conventional ALD and being compatible with roll-to-roll and other low cost and high throughput processing methods. Such an approach presents a method to enhance stability of flexible and transparent electrodes in a way compatible with industrial needs. The stability improvement observed is directly related to the thickness of the ZnO coating: the thicker the ZnO coating layer the better the thermal and electrical stability of the network.

However, ZnO coated AgNW networks are slightly less transparent when the ZnO coating thickness is increased, and thus a compromise in terms of ZnO thickness should be determined depending on the targeted application. The coated AgNW networks were evaluated as transparent heaters, proving to be much more stable than bare networks and thus reaching temperatures as high as 350 °C during a voltage ramp. In particular, the electrical failure voltage (and associated temperature) increased from 8 V (150 °C) to 14 V (350 °C) for bare and coated AgNWs, respectively, which represents a 133% increase in achievable temperature. Thanks to a simple model we showed that the origin of the stability enhancement can be associated with hindered atomic diffusion on the surface of the nanowires, since it is partially delayed due to the low diffusion of Ag atoms through the ZnO coating, thus preventing junction deterioration. The model enables to understand the evolution of network stability with coating thickness. Although the present work has been devoted to AgNW networks, the main physical ideas could as well be transferred easily to other nanomaterials.

47

 47

Figure 1 .

 1 Figure 1. (a) Scanning electron microscopy (SEM) image of the AgNW network deposited by spin coating and annealed at 250 o C for 30 minutes to sinter the junctions, (b) SEM image

  Figure 2a that for bare AgNW networks most of the nanowires are removed, showing the

Figure 2 .

 2 Figure 2. (a-b) Optical microscope images of the bare and coated (30 nm) AgNW network,

Figure 3 .

 3 Figure 3. a) Relative variation of electrical resistance (R/R 0 ) of bare and ZnO coated AgNW networks measured during thermal ramps from room temperature to 500 o C with a

  shows the voltage failure dependence versus ZnO thickness. This clearly shows that the ZnO coating enhances efficiently the electrical stability of AgNW networks. As in the case of the thermal stability enhancement provided by the coating, the electrical stability enhancement does not show a linear relationship with ZnO thickness.

Figure 4 .

 4 Figure 4. (a) Variation of electrical resistance for bare and ZnO coated AgNW networks

  2 10 10 A/m 2 for AgNW coated with ZnO. Such values of electrical current density are associated at the network level where enough AgNW have undergone a failure. These values at the network level can be compared with the one observed by Stahlmecke et al., 35 10 10 A/m 2 , for which void formation within an individual

  Figure 5e exhibits a schematic representation of the coated nanowires studied here (Fig 1c and 1d) while Figure 5f shows a SEM image of a AgNW coated with 25 nm ZnO. Based on the above arguments dealing with atomic surface diffusion, a simple way to estimate the effects of the ZnO coating on the AgNW network stability can be considered in terms of failure time . Let us consider the electrical stability and the effect on it of a ZnO coating (of a certain thickness L ZnO ) as reported in Figure 4. The failure time will be directly associated to a failure voltage during a voltage ramp through:

Figure 5 .

 5 Figure 5. Schematic representation of a junction between two adjacent AgNWs for a) as-

Figure 6 .

 6 Figure 6. (a) Electrical resistance of AgNW networks measured versus time during voltage

Figure 7 .

 7 Figure 7. Measured temperature versus voltage during a voltage ramp of 0.1 V/min for bare

Acknowledgements

This work was supported by the European Joint Doctorate FunMat (H2020-MSCA-ITN-2014, Project ID 641640) and by Campus France for the stay of Afzal Khan. The authors would also like to acknowledge the European Community through the funding support of the H2020-ICT 29-2014 LEO project and as well French National Research Agency in the framework of the "Investissements d'avenir" program (ANR-15-IDEX-02) through the project Eco-SESA. DMR acknowledges funding through the Marie Curie Actions (FP7/2007-2013, Grant Agreement No. 631111). VHN thanks the "ARC Energy Auvergne-Rhône Alpes" for economic support for PhD grants. This work was performed within the framework of the Centre of Excellence of Multifunctional Architectured Materials "CEMAM" n° ANR-10-LABX-44-01, through the Earth project. NDN acknowledges financial support by the European Joint Doctorate FunMAT program (H2020-MSCA-ITN-2014, Project ID 641640). The authors would like to warmly thank the following colleagues for fruitful discussions: Y. Bréchet, C. Celle, M. Lagrange, D.P. Langley, H. Roussel, T. Sannicolo and J.P. Simonato, as well L. Rapenne and E. Bellet-Amalric for their contribution to the nanowires characterization by TEM and SEM, respectively.

Silver nanowires (AgNW) coated with conformal thin ZnO layers by atmospheric pressure spatial atomic layer deposition (AP-SALD). The coating causes a drastic enhancement of the thermal and electrical stability of the networks which constitute very stable transparent electrodes. A model predicts that such stability improvement is the result of hindered silver atomic diffusion due to the presence of the thin oxide layer.