Guillaume Pirot

Emanuel Huber

James Irving

Niklas Linde

Reduction of

published or not. The documents may come

Introduction

Reliable predictions of groundwater flow and contaminant transport require adequate characterization of subsurface properties and their connectivity (e.g., [START_REF] Gómez-Hernández | To be or not to be multi-Gaussian? A reflection on stochastic hydrogeology[END_REF][START_REF] Zinn | When good statistical models of aquifer heterogeneity go bad: A comparison of flow, dispersion, and mass transfer in connected and multivariate Gaussian hydraulic conductivity fields[END_REF]. In this regard, limited number of data and knowledge of the field site implies that multiple geological conceptual models must be initially considered. That is to say models with different geometrical characteristics of the deposits, such as channels, lenses or layers. A general approach to compare alternative geological conceptual models is to perform Bayesian model selection based on field data acquisition and inversion.

It aims at estimating the Bayes factors, that is, the ratios of the estimated evidences (i.e., the integral of the likelihood over the prior probability density function) for the considered scenarios [START_REF] Kass | Bayes factors[END_REF][START_REF] Schöniger | Model selection on solid ground: Rigorous comparison of nine ways to evaluate Bayesian model evidence[END_REF]. However, reliable evidence estimators are costly because they necessitate a very large number of numerical evaluations of property models. As a result, modelers often assume a single conceptual model [START_REF] Ferré | Revisiting the Relationship Between Data, Models, and Decision-Making[END_REF] on which they perform inversion on the distribution of physical properties such as hydraulic conductivity, porosity or storativity [START_REF] Carrera | Estimation of aquifer parameters under transient and steady state conditions: 1. Maximum likelihood method incorporating prior information[END_REF][START_REF] Højberg | Model uncertainty-parameter uncertainty versus conceptual models[END_REF][START_REF] Eaton | On the importance of geological heterogeneity for flow simulation[END_REF] for a given geological conceptual model. The main risks associated with such a practice is underestimation of uncertainty and biased parameter distributions and predictions. There is, thus, a need for efficient, albeit more approximate, ways to compare alternative geological scenarios without resorting to formal evidence computations.

To enable comparison of geological conceptual models using a reduced number of costly forward simulations, [START_REF] Park | History matching and uncertainty quantification of facies models with multiple geological interpretations[END_REF] draw property models from each of the considered scenarios and calculate their data response. They then use multi-dimensional scaling (MDS) to reduce dimensionality, followed by adaptive kernel smoothing to estimate the probability of each scenario by comparing its distance to the reference data. Sometimes, it can be beneficial to base such comparisons on data types other than classical hydrogeological data (Huber and Huggenberger, 2016). Non-invasive geophysical data, for example, can provide substantial information about connectivity, structure dimensions and orientations, and thus might help to reduce geological conceptual model uncertainty. Notably, geophysical images reflect the sensitivity of the employed method to subsurface property variations. Thus, they can provide information about length scales and orientation characteristics of significant property boundaries. The wide range of available geophysical techniques offer flexibility to adjust resolution or depth of investigation, and to maximize the sensitivity to subsurface properties of interest [START_REF] Hubbard | Introduction to Hydrogeophysics[END_REF]. For instance, comparisons of seismic images [START_REF] Scheidt | Probabilistic falsification of prior geologic uncertainty with seismic amplitude data: Application to a turbidite reservoir case[END_REF] or of electric resistivity tomography (ERT) images [START_REF] Hermans | Uncertainty in training image-based inversion of hydraulic head data constrained to ERT data: Workflow and case study[END_REF] offer possibilities to falsify scenarios or reduce conceptual model uncertainty.

Possibly the simplest way to quantitatively compare geophysical images is to use a distance based on pixelwise (one-to-one) local comparison [START_REF] Hermans | Uncertainty in training image-based inversion of hydraulic head data constrained to ERT data: Workflow and case study[END_REF]. However, by using a local comparison, the probability of sharing a majority of similar pixel values and thus to observe small distances is quite low. So, when the main interest lies in the comparison of patterns and not the specific locations of property values, approaches relying on global geometrical characteristics are better suited. Approaches to sort and classify images in this way has been widely studied in the field of image processing [START_REF] Smeulders | Content-based image retrieval at the end of the early years[END_REF]. Among many alternatives, those based on discrete wavelet transforms have proven efficient to identify the images that are the closest in a large database. [START_REF] Suzuki | A distance-based prior model parameterization for constraining solutions of spatial inverse problems[END_REF] and [START_REF] Scheidt | Representing spatial uncertainty using distances and kernels[END_REF] use a distance based on wavelet decomposition [START_REF] Mallat | A theory for multiresolution signal decomposition: the wavelet representation[END_REF] of geological realizations for different scenarios to represent spatial uncertainty. [START_REF] Scheidt | Probabilistic falsification of prior geologic uncertainty with seismic amplitude data: Application to a turbidite reservoir case[END_REF] further apply this type of metric on seismic images to update probabilities of alternative prior scenarios. Nevertheless, distances based on wavelet decomposition rely on the comparison of coefficient histograms, which might hide spatial characteristics such as pattern connectivity.

It is, thus, important to also consider other distances, for instance, based on multiple-point histogram [START_REF] Boisvert | Multiple point metrics to assess categorical variable models[END_REF] or connectivity analysis [START_REF] Allard | Connectivity metrics for subsurface flow and transport[END_REF][START_REF] Meerschman | A practical guide to performing multiple-point statistical simulations with the direct sampling algorithm[END_REF], that allow quantitative comparison of the global spatial characteristics of interest obtained from field data with those obtained from synthetic modeling based on various scenarios.

So far and to the best of our knowledge, quantitative approaches to reduce conceptual geological model uncertainty using image comparisons did not consider multiple distance types and there has been no such application to GPR data. Traditionally, GPR data are interpreted qualitatively and its quantitative integration in subsurface modeling is largely unexplored. In the continuity of previous related works [START_REF] Park | History matching and uncertainty quantification of facies models with multiple geological interpretations[END_REF]Pirot et al., 2014;[START_REF] Scheidt | Probabilistic falsification of prior geologic uncertainty with seismic amplitude data: Application to a turbidite reservoir case[END_REF][START_REF] Hermans | Uncertainty in training image-based inversion of hydraulic head data constrained to ERT data: Workflow and case study[END_REF], we propose to extend such approaches to GPR reflection sections, using different distance measures of global geometrical characteristics. The three types of distances considered herein for the comparison of GPR reflection sections are based on 1) wavelet decomposition, 2) multiple-point histogram and 3) connectivity functions. In addition, the computed distances are analyzed and interpreted with a simple intuitive approach and with a more complex formal approach based on dimensionality reduction and mapping techniques.

The objectives of this work are i) to demonstrate how a simple but robust method enables the comparison of global characteristics of GPR reflection sections obtained from field data processing with those obtained from GPR reflection sections simulated from different scenario realizations; ii) to verify that GPR reflection sections can be used to reduce geological conceptual model uncertainty; iii) to investigate the relative strengths of three different distance measures for GPR data; and iv) to present follow-up strategies depending on the closeness or remoteness of simulated sections with reference sections obtained from field measurements.

To illustrate the proposed method, we consider GPR profiles acquired on the riverbed of the Tagliamento River, Northeast Italy [START_REF] Huber | Incorporating sedimentological observations, hydrogeophysics and conceptual knowledge to constrain 3D numerical heterogeneity models of alluvial systems[END_REF]. We consider three different geological conceptual models; each one of them being sub-divided in three sets of parameters (scenar-ios). For each of the nine resulting scenarios, 20 stochastic aquifer realizations are used as inputs for GPR simulations. The distances are used to produce a first ranking and to falsify unlikely scenarios. A dimension reduction technique called multi-dimensional scaling (MDS) followed by kernel smoothing are then used to estimate scenario probabilities.

The paper is organized as follows. Section 2 describes the distance measures considered and how they can be used to update scenario probabilities. Section 3 presents a fielddemonstration using GPR sections simulated from realizations of different geological conceptual geological models of the Tagliamento site (subsection 3.1). This section continues with the presentation of the migrated field GPR data and its processing steps (subsection 3.2), and ends with the simulation of migrated GPR profiles (subsection 3.3). Section 4 displays the results, which are further discussed in Section 5. Conclusions are given in Section 6.

Distances between geophysical images and estimation of scenario probabilities

In this section, we briefly review three distance measures that can be used to compare global geometrical characteristics of geophysical images. We then describe how approximate scenario probabilities can be obtained from field and simulated data through MDS and adaptive kernel smoothing [START_REF] Park | History matching and uncertainty quantification of facies models with multiple geological interpretations[END_REF].

Wavelet decomposition

One way to extract global characteristics of an image is wavelet decomposition [START_REF] Mallat | A theory for multiresolution signal decomposition: the wavelet representation[END_REF].

We consider in our work the same decomposition as [START_REF] Scheidt | Probabilistic falsification of prior geologic uncertainty with seismic amplitude data: Application to a turbidite reservoir case[END_REF]. Two geophysical images i 1 and i 2 are decomposed in two levels by a "Haar" wavelet [START_REF] Haar | Zur theorie der orthogonalen funktionensysteme[END_REF], which produces a series of coefficients (horizontal, vertical, diagonal and approximation) for each level. At each level, the histogram of each coefficient is discretized into bins b ∈ 1 . . . B, using the same binning for both images. For each level m ∈ 1 . . . M and each coefficient c ∈ 1 . . . C, a distance d JS between the two images is computed based on the Jensen-Shannon divergence between the probability distributions P m,c 1 and P m,c 2 derived from these histograms:

d JS (i 1 , i 2 , m, c) = d KL (P m,c 1 , P m,c 1 +P m,c 2 2 ) + d KL (P m,c 2 , P m,c 1 +P m,c 2 2 ) 2 , ( 1
)
where d KL (P, Q) is the Kullback-Leibler divergence between discrete probability distributions [START_REF] Kullback | On information and sufficiency[END_REF]. Then, the corresponding wavelet-based distance D w (i 1 , i 2 ) is:

P and Q computed as d KL (P, Q) = B b=1 P (b) log Q(b) P (b)
D w (i 1 , i 2 ) = M m=1 C c=1 d JS (i 1 , i 2 , m, c) M × C . (2)

Multiple-point histogram

Another way to quantify global spatial characteristics of an image is to define a summary statistic describing its multiple-point histogram [START_REF] Boisvert | Multiple point metrics to assess categorical variable models[END_REF]

D mph (i 1 , i 2 ) = M m=1 O o=1 |f o,m i 1 -f o,m i 2 | × (f o,m i 1 + f o,m i 2 ) 2 × M × O . ( 3
)

Connectivity measure

The final measure that we consider to quantify global characteristics of an image is connectivity [START_REF] Allard | Connectivity metrics for subsurface flow and transport[END_REF]. Indeed, subsurface property connectivity dictates subsurface flow paths and transport. Here we consider categorical geophysical images, but note, that the measure can be adapted to deal with continuous geophysical images (Pirot et al., 2014). We consider connectivity as the probability that two pixels belonging to the same class (a range of values) are connected, as a function of the distance and direction, similarly to the definition of a directional semi-variogram [START_REF] Matheron | Principles of geostatistics[END_REF]. By denoting C(i, a, l) the connectivity measure of a discrete image i along axis a ∈ 1 . . . A for a distance lag l ∈ 1 . . . L, the connectivity distance D c (i 1 , i 2 ) between discrete images i 1 and i 2 can be computed [START_REF] Meerschman | A practical guide to performing multiple-point statistical simulations with the direct sampling algorithm[END_REF] as

D c (i 1 , i 2 ) = A a=1 L l=1 |C(i 1 , a, l) -C(i 2 , a, l)| A × L . (4)

Estimation of scenario probabilities

To assess the probability of a scenario given a geophysical section, we follow the approach by [START_REF] Park | History matching and uncertainty quantification of facies models with multiple geological interpretations[END_REF]. Given a distance metric D and an ensemble of I images i, the distance between all pairs i j , i k of images define a dissimilarity matrix δ jk = D(i j , i k ). Multidimensional scaling (MDS, [START_REF] Cox | Multidimensional scaling[END_REF] is a method to represent the images as points in a low dimensional space, usually Euclidean. While principal component analysis (PCA) requires point coordinates, MDS can be used on data for which only the relative distances are known. This lower dimensional space is searched, such that the distances d jk between the points are as close as possible to the original dissimilarity matrix δ jk . MDS allows to map images in space, as points, for instance in 2D if using the two main dimensions. Now, we consider reference points related to reference images and a cloud of points related to images derived from a scenario. We can approximate the density of the cloud at any location of the low dimensional space, using adaptive kernel smoothing [START_REF] Ebeling | ASMOOTH: a simple and efficient algorithm for adaptive kernel smoothing of two-dimensional imaging data[END_REF]. For each scenario s, the density at one or several reference points (in the low dimensional space) can be computed as a scalar ρ s . The updated probability P of scenario s can then be approximated as P (s) = ρ(s) s ρ(s) . These updated probabilities are relative to the ensemble of considered scenarios, with P (s) the probability that an image generated from scenario s is the closest to the reference image.

Field application and GPR modeling

A pre-requisite to compare field and simulated data (Figure 1) is to apply equivalent data processing [START_REF] Hermans | Uncertainty in training image-based inversion of hydraulic head data constrained to ERT data: Workflow and case study[END_REF], but this is rarely sufficient because actual field conditions always differ from numerical implementations. Indeed, results obtained from the processing of geophysical data are prone to errors (e.g., [START_REF] Linde | Falsification and corroboration of conceptual hydrological models using geophysical data[END_REF] related to field data acquisition, simplifications in physical modeling or consequences of numerical modeling such as numerical and geometrical approximations. For instance, seismic or GPR geophysical images obtained from field data might include false discontinuities and their interpretation in terms of continuous connected structures or interface delineation necessitates expert knowledge. On the contrary, seismic or GPR geophysical images obtained from forward modeling, might reproduce property (dis)continuities too well and appear too clean to be representative of what would be expected for real data. To further reduce the remaining gaps between the results obtained from field data and from synthetic scenarios, it is necessary to include fit for purpose filtering [START_REF] Green | A transformation for ordering multispectral data in terms of image quality with implications for noise removal[END_REF][START_REF] Panagiotakis | Curvilinear structure enhancement and detection in geophysical images[END_REF] such that geophysical sections are not dominated by details/aspects that we do not seek to reproduce.

Study site and geological conceptual models

The study site considered is a portion of a sandy-gravel aquifer located near the city of Flagogna, Italy, within a portion of the active bed of the gravelly braided Tagliamento river (Figure 2). The Tagliamento river flows in the Friuli Venezia Giulia region, northeastern Italy, from the Carnian Alps to the Adriatic Sea. As the Tagliamento river is one of the few remaining large semi-natural rivers in the Alps [START_REF] Ward | A reference river system for the alps: the fiume tagliamento[END_REF] it was chosen as a study site to characterize the link between the topography of the active river bed and subsurface properties [START_REF] Huber | Morphological perspective on the sedimentary characteristics of a coarse, braided reach: Tagliamento River (NE Italy)[END_REF]. GPR data acquisitions and interpretations allowed to improve the characterization of scours and to model them (Huber et al., 2016).

In addition to improving the understanding of deposition and erosion processes (Huber and Huggenberger, 2016), this work inspired modelers to develop new methods, such as a pseudogenetic approach to produce heterogeneous models of braided-river aquifers (Pirot et al., 2015).

Assuming a braided-river type of aquifer, we wish to investigate which geological conceptual model is best suited to represent the porosity field. To this end, we consider subsets of reflection GPR sections in the saturated zone. Indeed, below the water table, GPR responses are strongly dependent on the porosity variations in the subsurface [START_REF] Daniels | Ground penetrating radar[END_REF]. We consider three different types of conceptual models of porosity , similar to those considered by Pirot et al. (2015) in their assessment of the impact of geological conceptual models on 3). Here we further assume that the braided-river aquifer is composed of three structural elements: gray gravel (GG), bimodal (BM) and open-framework (OW) deposits. Each distinguishable geobody or sedimentary deposit is a assigned a randomly drawn value from the porosity distribution, related to its structural element [START_REF] Jussel | Transport modeling in heterogeneous aquifers: 1. statistical description and numerical generation of gravel deposits[END_REF], as described in Table 1. The models are characterized by a horizontal discretization of 0.25 m and a vertical discretization of 0.01 m. The first geological geological conceptual model is represented by realizations from a pseudo-genetic (PG) algorithm (Pirot et al., 2015), which mimics deposition and erosion steps by stacking successive simulated topographies, and by imitating sandy-gravel material transport and deposition. Here, the main layers are populated with GG elements and the resulting cross-stratified deposits by successive BM and OW elements. A second geological geological conceptual model is a truncated multi-Gaussian (MG) model [START_REF] Emery | Tbsim: A computer program for conditional simulation of three-dimensional gaussian random fields via the turning bands method[END_REF], in which the locations above the highest threshold are populated with OW elements, the location between the two thresholds are defined as BM elements, and the remaining matrix is populated with GG elements. The third geological geological conceptual model is an object-based (OB) model (Huber et al., 2016) mathematically defined as a compound marked Strauss process. The OB simulates the formation of spoon-shaped structures on the river bed and the subsequent deposition of sediments over the whole river bed. The spoon-shaped structures are modeled by truncated ellipsoids with an internal OW-BM crossbedding and the sediments deposited on the river bed by horizontal layers of GG (e.g., [START_REF] Beres | Using two-and threedimensional georadar methods to characterize glaciofluvial architecture[END_REF][START_REF] Huggenberger | A sedimentological model to characterize braided river deposits for hydrogeological applications[END_REF]. The parameters underlying each scenario are summarized in Table 2; they were chosen to approximate the dimensions of scours that were estimated from field observations and from interpretations of migrated GPR sections.

GPR data acquisition and processing

The reflections in the processed and migrated GPR sections provide indirect information about characteristic geometric features. Such sections are used herein to compare, based on various global distance measures, different types of geological conceptual models. Five GPR profiles (REF01 to REF05) were acquired on the Tagliamento riverbed, orthogonally to the main flow direction. REF01 section is used for comparison with simulated data, while REF02 to REF05 are used to assess on-site data variability. The GPR data were acquired with a PulseEkko Pro GPR system (Sensors & Software Inc., Mississauga, Canada) using 100 MHz antennas and a measurement spacing of 0.25 m. A common mid-point (CMP) was performed to estimate the mean GPR velocity. The data processing steps are described in Table 3 and they were carried out with the RGPR package [START_REF] Huber | RGPR: a free and open-source software package to process and visualise ground-penetrating radar (GPR) data[END_REF]. The migrated section corresponding to the REF01 profile is presented in Figure 3a.

The processed migrated sections are thresholded into binary images to focus on the predominant aspects of the reflections. The amplitude of the processed GPR reflection section is similar throughout the image after applying the automatic gain control. Consequently, at all interfaces where porosity changes, the signal amplitude is similar, independently of the porosity contrast. We consider the first (negative) and last (positive) quartiles of the signal amplitude in the section. We retain the last quartile of the reflections (positive amplitude)

to define Class 1. Tests (not shown) indicated that it was not necessary to retain the first quartile (negative amplitude) to define another class, as the corresponding class would have almost the same geometrical characteristics as those of Class 1. Therefore, we use amplitudes below the 75 th percentile to define Class 2 (Figure 5a).

Table 3: Processing steps applied to field GPR reflection data.

Step

From aquifer porosity models to GPR reflection sections

In order to estimate the distances of each scenario realization to the reference GPR sections REF01, GPR reflection sections are simulated from the corresponding 2D porosity sections.

The processing steps are:

1. Realization of a facies/porosity model according to a geological conceptual model (scenario) as described in Section 3.1.

2. Porosity fields are converted into electrical property fields and velocity fields using the model by [START_REF] Pride | Governing equations for the coupled electromagnetics and acoustics of porous media[END_REF]. The petrophysical parameters (cementation index m, and dielectric constant of solid grains κ s ) are calibrated, such that the mean velocity of the corresponding porosity field is the same as the one used for the field data migration (vel = 100 m/µs).

3. Construction of a perfectly migrated GPR section (following the method developed by [START_REF] Irving | Inversion for the stochastic structure of subsurface velocity heterogeneity from surface-based geophysical reflection images[END_REF] by convolution of the propagated wavelet with a Primary Reflectivity Section. The propagated wavelet is estimated from field data processing step 5

(according to the method by [START_REF] Schmelzbach | Efficient deconvolution of ground-penetrating radar data[END_REF]. The Primary Reflectivity Section is derived from the previously obtained velocity model. A simple Gaussian horizontal filter is applied on the convolution result, to account for the Fresnel zone and whose width is determined by the dominant signal wavelength.

4. To mimic the effect of a constant velocity migration, the GPR reflection section generated with the actual velocities predicted from a porosity model is converted in the time domain before being back transformed into the depth domain using the same mean velocity as the one used in the migration of the field data (vel = 100 m/µs), and finally re-interpolated over a regular grid on the vertical axis.

5. 1D vertical Gaussian filter to slightly smooth the propagated wavelet with the same parameter as the one applied in the processing of the field data.

6. Automatic gain control to balance signal amplitudes with the same parameters as the one applied in the processing of the field data.

The resulting synthetic GPR sections (Figures 4b-j) are thresholded into binary images in the same way as the field data. The binary images resulting from the porosity images in Figures 3b-j are given in Figures 5b-j.

Results

For each of the three types of geological conceptual models and each of the three corresponding parameter sets (i.e., the nine considered scenarios), we generated 20 porosity realizations.

This means, that a total of 180 binary images were available for comparison with the binary reference section REF01 (Figure 5a). Wavelet-based, multiple-point histogram, and connectivity distance measures were computed between all possible pairs of field and synthetic binary images as follows. The wavelet-based distance uses B = 50 bins and M = 2 The distances obtained between all binary images and the Tagliamento reference section REF01 are displayed and grouped for each distance type by geological scenario (Figure 6).

To indicate the internal variability of the distances between the actual field data, the distances between binary reference section REF01 and other binary reference sections (REF02 to REF05) are gathered in a group denoted "REF". An acceptance threshold is defined by multiplying by 1.2 the maximum REF distance value. This subjectively-chosen acceptance threshold is used to select realizations whose distances to REF01 is similar to those of the reference sections. To better understand the generally-better performance of the PG-family as judged by the connectivity-based distance, we present connectivity functions in Figure 7 for some of the sections displayed in Figure 5. For the Class 1 components, the horizontal connectivity function (Figure 7a) is best reproduced by PG2, while the connectivity is overestimated for MG2 (by ≈ 0.08) and severely overestimated for OB1 (by 0.1 to 0.5). The vertical connectivity function (Fig. 7b) is best reproduced by MG2, while it is slightly too high for PG2 (at most ≈ 0.05 between 0.4 m and 0.9 m) and far too high for OB1 (by 0.1 to 0.2). For the horizontal and vertical connectivity functions of the Class 2 components (Figure 7c-d),

MG2 is found to reproduce them the best, while the connectivity is slightly lower for PG2

(by ≈ -0.02) and much too small for the OB1 scenario (up to -0.2).

To highlight the relationships between the distance types, we display three scatter plots (Figure 8). For each distance measure considered, the distances for all pairs of images are used to estimate the density of each scenario in the low dimensional space obtained by MDS. To estimate the updated probability of each scenario (Table 4), we limit the number of dimensions used such that 95% of the information is recovered. To achieve this, the two first MDS dimensions are sufficient for the wavelet-based distance, 14 are necessary for the multiple-point-based distance, and three are enough for the connectivity-based distance. For each distance, the estimated probability for a given scenario is proportional to the density of the cloud composed by the scenario realizations at the location of the reference section REF01 in the MDS space.

It informs about the probability that a realization from a scenario is closer to the reference section REF01 relative to the considered scenarios. Considering the wavelet-based distance, with an estimated probability of 85.9%, PG2 is the most probable scenario and MG1 is the second most likely one (14.1%). For the multiple-point histogram distance, PG2 is by far the most probable scenario (99.9%). For the case of the connectivity-based distance, MG2

is judged more likely (47.9%) than PG2 (33.6%) followed by PG1 (12.7 %), because it has fewer high and also the smallest distance value. If we average the probabilities over the types of distances considered, the scenarios that produce realizations that are the closest to the Tagliamento reference section REF01 is PG2, followed by MG2.

Discussion

Geological scenario falsification at the Tagliamento study site

By using three different distance metrics quantifying the agreement between field and simulated GPR sections, we reduce geological conceptual model uncertainty at the Tagliamento site. The direct analysis of the distances (Figures 6 and8) and the estimated probabilities for each type of distance (Table 4) led to similar conclusions. For the nine scenarios considered, two are judged significantly more suitable than the others: the PG2 scenario is the most suitable (its realizations are the closest to the Tagliamento reference section REF01), followed by the MG2 scenario. For both the wavelet-based and multiple-point histogram distances, PG2 is the most probable scenario. In the case of the connectivity-based distance, MG2 is judged the most probable scenario, followed by PG2.

To understand these rankings, let us consider the binary reference section (Figure 5a).

It reveals that: i) Class 1 components (main reflectors) have very small, small, intermediate and long length scales; ii) Class 1 components are sub-horizontal, and smaller components might present a stronger dip; iii) the interface between Class 1 and Class 2 components is irregular; iv) Class 2 components form a connected matrix. For the PG2 scenario realizations (Figure 5c), characteristic (i), (ii) & (iv) are present, but the interfaces are smooth. For MG2 scenario realizations (Figure 5f), characteristic (i), (iii) & (iv) are present, but the Class 1 components are too horizontal. The fact that scenarios PG2 and MG2 realizations fulfill three of these four visual criteria might explain the acceptable distance of their realizations to the Tagliamento reference section REF01. For the OB3 scenario realizations, none of the four criteria is fulfilled, which results in high values for all types of distance measures. From these results, it seems that the representation of cross-stratified deposits, interface roughness, and partially disconnected interfaces are important to reproduce reflection GPR sections at the Tagliamento site.

None of the proposed OB scenarios match the Tagliamento reference section REF01. We see two main possible explanations: 1) the geometrical parameters of this conceptual model were not well chosen, that is, the size of the scours and the layer thickness might be too large, the density of scours too small, the inner structure of the scours (i.e. inside the truncated semi-ellipsoids) have too thick deposits, when compared to the PG scenarios; or 2) this conceptual model is inherently unsuitable for this site (e.g., interfaces at porosity changes are too clean, without any contour irregularities or apparent roughness when compared to MG scenarios). This discussion also highlights that identifying the main characteristics present in the reference images and analyzing their absence or presence in images derived from various scenarios may help to propose new conceptual models or scenarios. This suggests a possible iterative process in which initial results are used to guide improvements in the conceptual models considered.

Comparison of distance measures

We now interpret our results to identify which distance-types are the most suitable. We observe a piecewise linear relationship between the wavelet-based distance and MPH-based distance (Figure 8a). Since there is less overlap between scenarios along the wavelet-based distance axis (Figure 8a-b), we conclude that it is more suitable than the MPH-distance to rank geological conceptual models and, to a lesser extent, their most appropriate parameters. However, the MPH-based distance is also able to classify models according to their geological conceptual model and scenarios (Figure 8a andc), but it performs less well than the wavelet-based distance to distinguish scenario PG3 from OB scenarios. This distance appears to better account for local structures (similar patterns between PG and OB) while the wavelet-based distance better accounts for global structures (different shapes: truncated ellipsoids versus the structures of PG models). Indeed, MPS algorithms often have difficulties in reproducing large scale connectivity even when using multi-grid levels [START_REF] Strebelle | Conditional simulation of complex geological structures using multiplepoint statistics[END_REF][START_REF] Mariethoz | The Direct Sampling method to perform multiple-point geostatistical simulations[END_REF][START_REF] Rongier | Reproduction assessment of connected geobodies in multiple-point simulation[END_REF].

The connectivity-based distances differ the most from the other distances and they dis- is that the background Class 2 is less connected (see Figure 7).

Overall, the results suggest that the wavelet-based distance provides the best ability for scenario differentiation. The connectivity-based distance is also interesting because it adds information and helps refining the scenario selection. Moreover, the connectivity-based distance is particularly interesting if the final application includes transport simulations, whose outcome strongly depends on property connectivity. We also would like to point to previous work (Pirot et al., 2014), which showed that the MPH-based distance is more sensitive to the sign of property contrasts while wavelet-based distance is more sensitive to the magnitude of property contrasts. Other fit-for-purpose distances could be considered and global integrative distances, i.e. that combine multiple global distance types could also be useful.

Influence of ranking method and parameter choices

We have seen that scenario falsification can be performed either by direct analysis of the distances or by estimation of updated probabilities per scenario using MDS followed by adaptive kernel smoothing. On the one hand, direct analysis of distances requires several reference images to define a reasonable acceptance threshold. On the other hand, the estimations of updated probabilities per scenario necessitate the computation of distances for all pairs of images within the ensemble composed of reference and simulated images. Since this cost increases as the square of the number of images, this can become computationally very demanding. Furthermore, rankings and falsifications based directly on distances of scenario probability estimations are relative to the ensemble of considered scenarios. In addition, small distances do not imply that the scenario sections are "surrounding" or "containing" the reference section in a space mapping the sections (see Figure 9).

Each type of distance requires specific parameter choices. Wavelet-based distances are parameterized by the type of wavelet used (Haar in our case), by the number of decomposition levels (two here) and by the number of bins (50 here). We tested (not shown) the sensitivity to different wavelets (e.g., Daubechies, Coiflets, Symlets, Mexican Hat) and obtained similar (30). A number of three [START_REF] Zhang | Filter-based classification of training image patterns for spatial simulation[END_REF][START_REF] Straubhaar | An improved parallel multiple-point algorithm using a list approach[END_REF][START_REF] Straubhaar | Parallel multiple-point statistics algorithm based on list and tree structures[END_REF] or four [START_REF] Strebelle | Conditional simulation of complex geological structures using multiplepoint statistics[END_REF][START_REF] Dell'arciprete | Comparison of three geostatistical methods for hydrofacies simulation: a test on alluvial sediments[END_REF] we did not vary these parameters, because the connectivity functions (Figure 7) appears to be well defined.

Perspectives

In the presented case-study, we threshold the reflection GPR sections as part of the data processing (to focus on the main aspects of the reflectors) and limit our comparison to binary geophysical images. One could also apply the proposed methodology to continuous images.

It would then be straightforward to compute a distance based on wavelet decomposition.

However, multiple-point histograms and connectivity functions as defined in Section 2 are applicable to discrete domain images only. One solution is to threshold the continuous images, as we did here, in a reasonable number of classes, to retrieve and compare the most important features from the images. Of course, this implies some qualitative assessment of which features are the most important ones, depending on the target of modeling. Another possibility is to adapt the definition of the multiple-point histogram and of the connectivity functions, such that they can be applied to continuous images. For instance, we could rely on the definition of distances between continuous patterns [START_REF] Mariethoz | The Direct Sampling method to perform multiple-point geostatistical simulations[END_REF] and on the identification of pattern clusters to build a multiple-point histogram between continuous images; the pattern clusters could be referred to as the histogram bins, and a pattern could be assigned to the closest bin/cluster; it would though depend on the number of clusters and how they are identified. Regarding the connectivity-based distance, the simplest option would be to define the connectivity as a function of the threshold [START_REF] Meerschman | A practical guide to performing multiple-point statistical simulations with the direct sampling algorithm[END_REF][START_REF] Allard | Connectivity metrics for subsurface flow and transport[END_REF], as the probability that two pixels are both above or both below a threshold.

While migrated GPR sections obtained from field data are somehow affected by 3D geological heterogeneities, the simulation of GPR reflection sections is performed from 2D porosity sections and does not account for 3D effects. The binary thresholding is a way to focus on the reflections of interest and to reduce the impact caused by the inherent limitations of the forward modeling, such as considering 3D effects negligible, grid resolution, different coupling effect at the surface, non-horizontal antennae at all times due to small changes in topography, approximations of the propagated wave, estimation of the attenuation with depth, etc.

A consequence is that we loose some information about porosity contrasts. Here, it allows to simulate GPR reflection sections very efficiently, and thus to perform conceptual model uncertainty reduction. A way to account for 3D effects would be to perform full-waveform GPR modeling over 3D porosity models. It would tremendously increase the computational requirements, and consequently would make conceptual model selection and falsification very costly. However, characterizing the effects of such model simplifications could improve (quantitatively) our understanding of GPR modeling errors and allows us to mitigate potential bias effects.

Conclusions

We have demonstrated how global distances (defined from wavelet decomposition, multiplepoint histograms and connectivity analysis) between geophysical images allowed us to falsify seven out of nine considered geological scenarios at the Tagliamento site. By considering GPR sections from the Tagliamento aquifer, we find that cross-stratified deposits and irregular property interfaces are important features to reproduce. An underlying assumption of this work is that the results obtained by model comparison with geophysical data are informative for subsurface flow and transport. This assertion should be tested by tracer tests, that are, up to date, not available at the Tagliamento site. We have found that scenario falsification can be performed either by direct analysis of the distances or by estimation of updated probabilities. Direct analysis is faster, more intuitive and rely on the definition of a subjective acceptance threshold that is informed by the magnitude of distances computed between several reference sections. Computation of scenario probabilities using MDS to map geophysical images as points in a lower dimensional space, followed by adaptive kernel smoothing to estimate scenario probabilities, is more advanced and requires more computing resources. The use of distance comparisons in geophysics also serves to select new parameter sets or to propose new geological conceptual models, in order to further close the gap between simulated sections obtained from an initial set of scenarios and the reference sections. This approach can be used for any type of geophysical images, as long as the geophysical modeling and processing step can be simulated in an effective and trustworthy way. The most convenient distance of those considered is the wavelet-based distance, which is the fastest to compute and it offers the best clustering of scenarios. The connectivity-based distance add further independent information and should be considered if structure connectivity is expected to have an impact on the prediction variables of interest. This work proposes a way forward to use uninterpreted GPR data, in contrast to hand-drawn geological deposit interpretation, for quantitative subsurface characterization.

Figure 1 :

 1 Figure1: Overview of the workflow to reduce geological conceptual model uncertainty. On the left, the path of arrows represents field data processing; on the right, the three vertical arrow paths represents the workflow for three distinct scenarios; at the bottom, a red cross illustrates scenario falsification and a green mark indicates scenario compatibility.

Figure 2 :

 2 Figure 2: Site location in Italy (map from http://www.pedagogie.ac-aix-marseille.fr/ jcms/c_67064/en/cartotheque); position of the GPR profiles over an aerial photograph of the Tagliamento river, south east of Flagogna (Google maps satellite image).

Figure 3 :

 3 Figure 3: Example porosity sections for different geological scenarios that are to be compared to (a) a reference GPR reflection section processed from field data (REF01); (b), (c) & (d) porosity sections from pseudo-genetic model realizations for parameter sets PG1, PG2 & PG3, respectively; (e), (f) & (g) porosity sections from truncated multi-Gaussian model realizations for parameter sets MG1, MG2 & MG3, respectively; (h), (i) & (j) porosity sections from object-based model realizations for parameter sets OB1, OB2 & OB3, respectively.

  Figure 4: (a) Processed and migrated GPR reflection section from field data (REF01); (b), (c) & (d) GPR reflection sections simulated from pseudo-genetic porosity model realizations for parameter sets PG1, PG2 & PG3, respectively; (e), (f) & (g) GPR reflection sections simulated from truncated multi-Gaussian porosity model realizations for parameter sets MG1, MG2 & MG3, respectively; (h), (i) & (j) GPR reflection sections simulated from object-based porosity model realizations for parameter sets OB1, OB2 & OB3, respectively.

Figure 5 :

 5 Figure 5: Images obtained after thresholding the example sections represented in Figure 4. (a) The binary geophysical image obtained from field data; (b), (c) & (d) binary geophysical images obtained from pseudo-genetic porosity model realizations for parameter sets PG1, PG2 & PG3, respectively; (e), (f) & (g) binary geophysical images simulated from truncated multi-Gaussian porosity model realizations for parameter sets MG1, MG2 & MG3, respectively; (h), (i) & (j) binary geophysical images simulated from object-based porosity model realizations for parameter sets OB1, OB2 & OB3, respectively.

Figure 6 :

 6 Figure 6: Distance to Tagliamento reference section REF01; plots grouped by scenario for (a) wavelet-based distance, (b) MPH-based distance and (c) connectivity-based distance. REF denotes distances for other binary reference sections (REF02-REF05) with respect to REF01 and the red line corresponds to the acceptance threshold.

  A piecewise linear correlation between wavelet-based and multiple-point histogram distances is clearly visible in Figure 8a, in which a first segment corresponds to the PG and MG scenarios and a second to the OB scenarios. It also shows the ability of waveletbased and multiple-point histogram distances to distinguish between the different conceptual

Figure 7 :

 7 Figure 7: Example of connectivity functions for a selection of binary geophysical images (REF01, PG2 sim, MG2 sim & OB1 sim from Figure 5); (a) horizontal connectivity functions for Class 1 (white) components; (b) horizontal connectivity functions for Class 1 (white) components; (c) horizontal connectivity functions for Class 2 (gray) components; (d) horizontal connectivity functions for Class 2 (gray) components.

Figure 8 :

 8 Figure 8: Distance to Tagliamento reference section REF01 visualized as scatter plots grouped by scenario: (a) MPH-based distance as a function of wavelet-based distance; (b) connectivity-based distance as a function of wavelet-based distance; (c) connectivity-based distance as a function of MPH-based distance. REF denotes other binary images processed from additional GPR profiles (REF02-REF05) acquired at the study site and the red line corresponds to the acceptance threshold.

  play a weak log-linear relationship with the wavelet-based distances. They are effective in rejecting the MG1, MG3 and all OB scenarios. The connectivity-based distance clearly separate the OB models from the other model classes (as shown by Figures 6c and 8b-d) as the reflectors (Class 1) in the OB models are much too connected in length. A corollary of this

Figure 9 :

 9 Figure 9: Mapping of the simulated and reference sections in the two first dimensions of the MDS space; (a) for the wavelet-based distance; (b) for the multiple-point based distance; (c) for the connectivity-based distance.

  multigrid levels is commonly chosen to generate realizations with tree-or list-based MPS algorithms to capture patterns at multiple scales. The pattern geometry is a basic square which does not favor any anisotropy. The pattern size is kept relatively small to ensure the possibility to encounter similar patterns between images. A smaller pattern size (3 × 3 pixels window) was tested, but led to similar results. The number of most-frequent patterns is limited to 30 to avoid the comparison of single occurrences that are present only in one of the images. Increasing the number of most-frequent patterns would increase unnecessarily all distances. Decreasing the number of most-frequent patterns would reduce the distances between images. Connectivity-based distances are parameterized by investigated directions and lag width, similarly to the computations of semi-variograms. Here

Table 1 :

 1 Probability density function (pdf) properties of the porosity for each structural element (from[START_REF] Jussel | Transport modeling in heterogeneous aquifers: 1. statistical description and numerical generation of gravel deposits[END_REF].

	Structural Element Pdf Law Porosity Mean (%) Porosity Standard Deviation (%)
	GG	normal	20.1	1.4
	BM	normal	18.8	3.9
	OW	normal	34.9	1.4

Table 2 :

 2 Parameter choices for each scenario grouped by type of geological conceptual model.

	Scenario Example			Parameters
			Scalability	Scalability	Aggradation	Number of	Deposition
			Width	Depth	Range (m)	Iterations	Intensity
	PG1	Figure 3b	1	1	[0.05 ; 0.125]	8	5
	PG2	Figure 3c	1/2	1.6	[0.05 ; 0.125]	8	5
	PG3	Figure 3d	1/3	1	[0.2 ; 0.25]	6	3.5
			Variogram	Horizontal	Vertical	OW element BM element
			Model	Range (m)	Range (m)	Proportion	Proportion
	MG1	Figure 3e exponential	50	3	25%	25%
	MG2	Figure 3f exponential	25	0.5	25%	25%
	MG3	Figure 3g exponential	70	5	25%	25%
			Width	Width/height Layer Poisson Horizontal Strauss process
			Range (m)	Ratio	Process (λ)	β	γ
	OB1	Figure 3h	[10 ; 20]	[11 ; 18]	0.1	10 -3	0.5
	OB2	Figure 3i	[22 ; 33]	[11 ; 18]	0.1	5 10 -4	1
	OB3	Figure 3j	[35 ; 53]	[11 ; 18]	0.1	2.5 10 -4	1

Table 4 :

 4 Estimated scenario probabilities (%) computed for each type of distance by adaptive kernel smoothing on MDS representations of the simulated and reference sections; values smaller than 0.1% are not displayed; for each type of distance (row) the probabilities sums to 100%.

	Scenarios

Acknowledgments

The authors thank two anonymous reviewers for their constructive comments, and Céline Scheidt for sharing her MATLAB code to compute updated probabilities. The other MAT-LAB codes and the synthetic data used to test the proposed method are available upon request at guillaume.pirot@unil.ch (note that multiple-point histograms computations in that code require the Impala software, whose academic license is available upon request at philippe.renard@unine.ch). The software package to simulate the object-based model is freely available at github.com/emanuelhuber/CBRDM. The GPR data are available upon request at emanuel.huber@alumni.ethz.ch.


Valérie Gay
email: valerie.gay@uts.edu.au

Frédéric Weis
email: frederic.weis@irisa.fr

Creating smarter spaces to unleash the potential of health Apps

Keywords: Smart Space, Ubiquitous Applications, Edge/Fog computing, e-Health, Chronic disease management, Asthma Management, IoT, Context awareness

Technologies necessary for the development of pervasive health apps with intensive and seamless interactions with their environments are now widely available. Research studies and experimentations have demonstrated the real ability for health apps to interact with their environment. However, designing, testing and ensuring the maintenance and evolution of pervasive health apps remains very complex. In particular, there is a lack of tools to enable developers to design apps that can adapt to increasingly complex and changing environments (sensors added or removed, failures, mobility etc.). This paper reflects our vision to reduce this complexity and is based on our current research work on smart environment and personalized health monitoring apps. It uses SAM, a smart asthma monitoring app as an illustration to highlight the need for a comprehensive set of new interactions to help health app developers interact with the users' environment, and more specifically get a smarter access to the data. Some requirements can be on the minimum quality level of the data and the way to adapt to the diversity of the sources (data fusion / aggregation), on the network mechanisms used to collect the data (network/link level) and on the collection of the raw data (sensors). It discusses possible solutions to address these needs.

Introduction

Smart space eco-systems are becoming a mass market reality. Large areas of our living spaces are now instrumented: short/long-range and low energy communications, a broad variety of visible (smart objects) or invisible (sensors and actuators) objects. The applications of the smart spaces are boundless. There is a real ability to adapt the smart space to the behaviors and needs of users, however mainstream applications are barely existent.

In this paper we use the definition of [START_REF] Jayaraman | On-the-Fly Situation Composition within Smart Spaces[END_REF] for smart space: A smart space is a smart environment that has commuting embedded into it and can provide information that can be used to model the real world into the virtual computing world.

This paper envisions smarter spaces with massive but scalable interactions between individuals and their everyday working and living environments such as residential housing, public buildings and vehicles. It focuses on the needs of app developers for better ways to interact with their environment and in particular use the opportunity of tapping into the vast amount of data collected from the physical environments for health apps. Today, this can only be done if the environment has especially been developed for it (for example "smart" hospital room). They are difficult to adapt to increasingly complex situations, even though environments in which they evolve are more open, or change over time (new sensors added, failures, mobility etc.)

Being able to develop health apps that are ready to deploy and evolve in everyday environments should allow significant cost reduction. Unfortunately, designing, testing and ensuring the maintenance as well as the evolution of such a pervasive application remains very complex. In our view, the lack of resources by which data issued from real environments are made available to health app developers is a major concern.

This paper is looking into this lack of data exposure and its objective is to explore how application development, setup and maintenance could be simplified in an open smarter space eco-system. Using an existing health app we developped as a basis, section 2 presents an example on asthma management and shows how this chronic condition would benefit from access to the data issued by the smart space eco-systems. Section 3 gives our motivation for a new approach and focuses on one main challenge: providing health app developers with a better and smarter access to the environment data. The paper concludes with perspectives and future work.

Illustration: SAM -Smart Asthma monitoring

One chronic condition that would benefit from getting better data from smart space eco-systems is asthma. According to GAN 4 in 2014, asthma may affect more than 334 million people worldwide. A few projects focusing on Asthma management exist [START_REF] Votis | Analysis, modelling and sensing of both physiological and environmental factors for the customized and predictive self-management of asthma[END_REF] but they either use a tightly controlled environment or they do not adapt to their environment and rely mainly on personal area networks and specialized individual sensors (e.g. peakflows). Some cloud infrastructure solutions exist [START_REF] Huang | Allergy and asthma care in the mobile phone era[END_REF]. However they do not have the level of pervasiveness this paper envisions.

Figure 1 illustrates our vision. It features one citizen, Alex, an active middle age man that is prone to asthma attacks mostly triggered by pollen and pollution and a team of clinicians that help him managing his condition. Alex is using the app we developed. The app collects health and fitness data from different sources and gives Alex the possibility to aggregate data from different sources. The app interacts with a wide range of wireless devices and wearable health trackers and also aggregates data from third party apps. It offers Alex personalized exercise tracking and monitoring of biometric data such as heart rate, respiration, body temperature, weight, food intake, blood pressure, cholesterol, asthma, blood glucose and many more. Three important steps are involved in his asthma management in the smart space: Data collection, Data analysis and Nudge, feedback and service presentation.

Cloud Local environment

Current

Data collection for asthma

On a typical day, while Alex is living an active life, the smart space collects as much relevant data as possible in order to offer Alex the best overview on his asthmatic condition. His app collects data from wireless devices supporting open standard protocols such as the Google Android Wear smart watches and fitness trackers that allow third party developers to retrieve the data directly from the device (Figure 2, box 1). It also collects data from devices (peakflow meter, pollution meter) paired with Alex's mobile device (Figure 2, box 2) via Bluetooth or ANT+. Alex is in control of the data he wants to collect and how it is collected. He can enter data manually (Figure 2, box 3) and annotate a reading by adding comments and contextual information, such as extra performed activities. For closed and proprietary wireless devices (Figure 2, box 4 and 5), data is obtained via their server through an open API (e.g Fitbit, Jawbone, Withings, myFitnessPal and Fatsecret).

Data is also collected from the smart space in order to offer more fine-grained analysis (Figure 2, box 6). The context of Alex can be important and one of the challenges is to identify and collect the data that is important for a particular user in order to provide timely feedback. For Alex, our app captures locationbased information about the Air Quality Index (AQI) and pollen levels. In the global scenario summarized in Figure 1, our mobile application collects the data within Alex's smart space and in the cloud, for example, to get real-time data on the pollen and pollution levels in the different locations Alex intends to travel to today. Another example is to get forecast on air quality in Alex's smart home environment.

Fig. 2. Our app's eco system 2.2 Data analysis and data storage for asthma Currently, the app stores all data locally on the mobile device, or on the cloud on a server chosen by the user. Alex opts for that in order to get more intelligent feedback based on his long-term history and possibly by comparing his data with the data of other citizens with similar conditions. Some of his sensors can also contribute to crowd sensing and help other asthma sufferers. Alex also opts for his clinicians to get access to his data and decide, for example, that there is a need to make adjustments to his asthma plan. This option gives Alex's clinicians an overall view on his health and fitness data and help them making a more accurate and personalized assessment. This could also allow hospital to estimate how many hospital beds they need to reserve for victims of asthma attack.

Nudge, feedback and service presentation

For the feedback, we opted for a user-centred approach where Alex is in control of the way the data is presented to him. Triggers are personalised to Alex and the app provides threshold to trigger an alert. For example, the app would alert Alex when the air quality deteriorates and gets under his thresholds and raises an alarm if it gets hazardous for his condition. A similar approach will be taken for the pollen level. Alex can either see the data in real-time, or just get a notification or nudge when necessary or even ignore the data and get a phone call from his doctor if he needs to do something. Alex opted for a nudge when needed and an automatic reaction from his smart space when possible (automatic closing/opening of doors and windows at home or in the car, adjustment of the air conditioning). He selected to have real-time analysis and get warnings and references to his asthma plan based on his personal threshold levels on several parameters (including pollen and pollution level, peak flow readings). He may also be contacted by his clinicians to make an appointment or to adjust his asthma plan. Based on his calendar that contains the place and activities, he is planning to go to, he may also receive some advice to change plans (e.g.: work from home, exercise indoor) based on the location's current and forecasted data. Alex lives in Australia and he has an asthma action plan: a written plan, prepared for him by a healthcare professional. It allows him to reduce the severity of acute asthma flare-ups. Asthma action plans help recognize worsening asthma and advise what to do in response. They typically include: a list of prescribed medication, including dosage; what to do when asthma gets worse; what to do in an emergency; the name of the health professional predominantly responsible for treating the person and the date. In the future we would like Alex's asthma plan to be integrated in a health applications that would monitor the different triggers in the smart spaces and react according to the patient's asthma plan.

Challenges to develop SAM to support this vision

Most of the vision presented in this example is feasible today but parts of it still bring interesting challenges. This section highlights some of the challenges developers have to be able make this vision possible.

Silos There is an increasing amount of health-and fitness-related information that has been collected and stored in the cloud. However, the data usually resides in silos and in addition, in most cases health and fitness data is separated. According to Mandi et al [START_REF] Mandl | Driving innovation in health systems through an apps-based information economy[END_REF], these data streams will initially remain confined to their respective platforms and will have very limited ways to integrate with electronic health records (EHRs). To provide better health outcomes and better patient engagement, a complete picture is needed which combines informal health and fitness data collected by the user, together with official health records collected by health professionals. Combining these two streams, the data can be analyzed using data analytics and health professional expertise to offer better personalized advice and care. There is good evidence that the integration can improve therapeutic management [START_REF] Becker | mhealth 2.0: Experiences, possibilities, and perspectives[END_REF].

For chronic disease patients, there is often a need to monitor several physiological parameters and not just their activity level or environmental data. The presence of data silos prevents these users from obtaining an overall view of their condition. To engage these users in active monitoring of their condition, it is important to have all health and fitness data in one place and give them a personalized overview and feedback on trends and progress. Moreover, device man-ufacturers develop their own protocol and data formats to retrieve data from the device or to send commands to the device. Some protocols are straightforward, using plain text to send or receive data. Many however, implemented complex protocols with numerous commands to control and exchange data. Without a detailed specification it is impossible to communicate with these devices.

Quality and privacy of the data With the increase of data sources comes the need to be able to differentiate the sources based on their quality, security and trust levels. There is a need to express requirement on the network mechanisms used to collect the data (network/link level) security aspects in particular -and on the collection of the raw data (sensors). This would lead to a better service for asthma sufferer as their asthma management would be able to dynamically adjust itself to the status of the smart-space (e.g. change the air conditioning at home or in a personal car).

Abstraction level Health app developers also need to get the right level of abstraction to trigger an action. To get to that level of abstraction (e.g. is the pollution level too high wrt the users thresholds) there is a need to combine data at different levels (e.g. raw or elaborated data).

Detection of unreliable sensors

The reliability of the sensors varies widely and reliability is very important for health apps like SAM and health professionals that need to make a diagnostic. In Alex's case, air pollution data comes from accurate and expensive sensors like the ones of his high tech building at work, but it also comes from the low cost gadget Alex is carrying. The gadget might be untested and therefore not be fit for purpose to justify, on its own, a change to his medication or his asthma plan. For example, devices made for the fitness market are not necessarily approved by the Food and Drug Administration (FDA) and as such even less reliable. If we are able to tags the source of the reading, it would be beneficial for a health professional in his/her assessment of the data quality. Many health professionals discard self-collected health. Fitness and environmental data due to the unreliability of the data and prefer to use their own data for diagnosis. We believe that over time more health professionals will accept the data if the source is properly tagged so that they know which device, or which application, generated the data.

Smarter access to the local data

The technologies supporting the pervasive services presented in the asthma management example, with intensive and seamless interactions between individuals and their environments, are now widely available. Spaces in which Alex lives are increasingly enriched with these technologies. Alex's asthma can be monitored in his working and living environment: residential housing, public buildings, vehicles, etc. as illustrated on Figure 3a. These environments are themselves organized in silos that do not interact with each other and do not coordinate well. In each environment, designing an application involves implementing one or more logical control loops which include four stages: (1) data collection in the real environment, (2) the (re)construction of information that is meaningful for the application and (3) for decision making, and finally, (4) action within the environment. While many decision-algorithms have been proposed, the collection and construction of a reliable and relevant perception of the environment still poses major challenges.

Developing Health apps that are ready to evolve in different environments and that overcome the existing silos would involve significant cost reduction. Unfortunately, as illustrated in section 2, designing, testing and ensuring the maintenance as well as the evolution of a pervasive health application remains very complex. In our view, the lack of resources by which properties of the real environment are made to application developers is a major issue. Most current solutions are based on a massive collection of raw data from the environment, stored and processed on remote servers [START_REF] Mun | Peir, the personal environmental impact report, as a platform for participatory sensing systems research[END_REF] and use specialized cloud services for delivering pervasive services. These solutions isolate the developer from the real environment that often needs to be depicted according to complex, heavy and expensive processes. However, new paradigms like edge or fog computing [START_REF] Shi | Edge computing: Vision and challenges[END_REF] offer alternative ways of bringing computing capabilities to the edge of the network. Objects and infrastructure integrated into user smart space often provide a more suitable support to pervasive applications. The objects can carry different kinds of data, like current smart space air quality. The description of the actual state of the environment can be richer, more accurate, and, meanwhile, easier to handle. The processing can be distributed by being built directly into the environment, facilitating scalability and resilience by the processing autonomy; and finally, moving the processing closer to the edge of the network avoids major concerns about data sovereignty and privacy encountered in cloud-based infrastructures. Figures 3b and1 illustrated a plausible "distribution" of an asthma management scenario between Alex's local environment and cloud-based services. The best place to use the data depends partly on what kind of data is needed for making decision and trigger processing. Fine decisions can be made closer to the objects producing the data. Local data characterization and local processing reduce the computing and storage resources in the cloud. The latter can be always used to store selected or transformed data for global historical analysis or optimization, for example for health and fitness history.

There is a need for a comprehensive set of new models and architectures to help smart space designers making those distribution choices. We propose a design approaches that exposes local properties of the environment to pervasive applications. The developer must be able to manage and to enrich locally data produced in the environment. In return the application would then be able to build their knowledge about their environment (perception) in order to adjust their behavior (eg. level of automation) to the actual situation. In the following, we present three research directions to develop models of data-centric service architecture.

Characterizing data quality

Several studies have proposed data quality analysis [START_REF] Jayawardene | An Analysis of Data Quality Dimensions[END_REF] and methodologies for identifying data quality issues [START_REF] Challa | A data driven approach for discovering data quality requirements[END_REF]. The quality of the raw data collected in the smart space (temperature, sound level, pollen level etc.) are of crucial importance to assess the quality of the perception of the environment and therefore to ensure correct behavior of a health pervasive app. Poor quality data is often a result of poor representation of the physical environment. Data quality can be characterized through several dimensions. Here the term dimension is defined as a "measurable extent of a particular kind": timeliness, accuracy, confidence, reliability, confidentiality, etc. For example, timeliness refers to time data should be available for use by an app within an acceptable time relative to its time of creation.

In an environment instrumented with smart objects and sensors interconnected with low energy communications, the way data is produced, consolidated, and aggregated has an impact on its quality dimensions. Moreover part of these quality dimensions requires gathering information from several communication layers (link and network layers) at several entities (the object itself, access point, router etc.) to be computed on the fly. Lightweight cross-layer interactions to collect relevant data have to be designed for this purpose, the goal is to expose "augmented" data (initial raw data + some pertinent dimensions selected by the app developer) to the pervasive app. Figure 4 illustrates these principles.

A frugality principle should guide the design, therefore it is not appropriate to build all imaginable attributes. It is necessary to identify attributes relevant to the application and to have mechanisms to activate or deactivate at run-time the process to collect them. A "frugal" service uses as few resources as possible, while having the possibility to use much more when it is required. This is quite a different approach than best-effort services where a networked application uses all the resources it can use and if there are not enough resources, it tries to mitigates the effects.

Raw data app

Exposition of data

Tools for data quality characterization

Physical env.

dimension Fig. 4. Enriching data collected in smart environment using cross-layer interactions

Providing pertinent abstractions

Raw data should only be used directly to determine low-level abstractions. Further help in abstracting from low-level details can be provided by data fusion.

The pervasive community increasingly understands that developing pervasive apps should be supported by such data fusion tools [START_REF] Talebifard | A data fusion approach to context-aware service delivery in heterogeneous network environment[END_REF]. While many definitions exist for data fusion, we use the one in [START_REF] Khaleghi | Multisensor data fusion: A review of the state-of-the-art[END_REF]: Information fusion is the study of efficient methods for automatically or semi-automatically transforming information from different sources and different points in time into a representation that provides effective support for human or automated decision making.

A good construction of a meaningful information (for example "air quality in the home", "level of activity in a room", "posture of a person" etc.) reduces the complexity and helps developers to concentrate on the application logic.

Moreover, the reactivity required in pervasive systems and the aggregation of large amounts of data (and its processing) in the cloud are antagonists. Fusion tools that can be deployed closer to the edge of the network are needed.

Traditional data fusion systems are often designed to use a set of dedicated sensors. Moreover, most approaches rely directly on sensor measures making them hard to deploy in many different environments and for different pervasive apps. They require an ad hoc tuning phase as the behavior of a sensor can change from one environment to another. Finally, the techniques to be used in our approach are guided by different criteria [START_REF] Pietropaoli | Multi-sensor data fusion within the belief functions framework -application to smart home services[END_REF]: (1) the relevance of abstractions produced for pervasive applications, (2) the anonymization of exploited raw data, (3) the processing time and computing power required for processing, (4) the ignorance of the system (the system should be able to be indecisive when it does not have enough clue of what is going on) and ( 5) the no predetermined sensor configuration (the data fusion process should be if possible independent from the configuration of sensors).

We proved [START_REF] Pietropaoli | Multi-sensor data fusion within the belief functions framework -application to smart home services[END_REF], through the use of dedicated algorithms and a layered architecture that data fusion tools can be deployed using the Belief Function Theory (BFT) [START_REF] Dempster | Upper and lower probabilities induced by a multivalued mapping[END_REF]. Raw data collected in the environment can been seen as an "accumulation of evidence". BFT enables the use of multiple sources to infer a meaningful abstraction defined by app developers. In our approach, an abstraction used by pervasive app is computed by a chain of three tasks as illustrated in Figure 5: -Model. The transition between a raw sensor value and a model of the sensor is made through the use of a "belief function" which maps each sensor value to a possible state (state 1, state 2 etc.) of the abstraction. A belief function represents the "degree of belief" associated to each possible value of the abstraction. -Fusion. Then a set of belief functions (corresponding to a set of sensors) can be combined (fused). -Decision. Finally the system can decide what is the "best" value for the abstraction.

Detecting new sensors / detecting unreliable sensors

One of the reasons why pervasive apps are not affordable for mainstream adoption comes from the necessity to control all the aspect of the environment to ensure the correctness of the decision-making. This is especially true for apps that move in various unknown and uncontrolled environment. Building a fit for purpose perception using available data from the environment is one way to simplify the deployment since it alleviates the need to know a priori the peers we need to communicate with. As stated in section 3.2, it is not possible to use data as it is produced; it is then necessary to calibrate the raw data issued from a sensor to take into account its specificities. Such a calibration depends on the usage, this is why it is often done upon a complex installation phase on the field and not during the manufacturing process. Some solutions exist to solve this issue. As an application could not know a priori all sensors, if nobody does it before, it has to establish knowledge about sensors and the quality of the information they generate (eg. how to calibrate them) for a set of purposes. Such knowledge could be stored in the environment, then the information could be provided to other applications or instances when needed.

More generally to ease the design of new applications and to align the development of new products with the ever faster standard developments, continuous integration could be used in parallel with continuous conformance and interoperability testing. New shared platforms like [START_REF]F-interop: Platform for online interoperability and performance test[END_REF] aim at facilitating this by providing remote testing tools. Unfortunately, from an application perspective it is not possible to be certain that all potential peers in the surrounding have a compliant behavior for the property it needs [START_REF] Calinescu | Self-adaptive software needs quantitative verification at runtime[END_REF]. Moreover, upon failure or security breach, a piece of equipment could stop to operate properly and lead to global misbehavior. Testing devices at runtime is then necessary. As for calibration the result of such conformance or interoperability tests could be stored safely in the environment by authoritative testing entities. Applications could then interact with the device with a higher confidence. The confidence level of a device could be part of the quality attribute of the information it contributed to generate.

Conclusion

We highlighted some of the complexity and the challenges of designing, testing and maintaining pervasive health app that need a great deal of interactions with smart spaces. An illustration based on asthma management highlights the importance of the interactions with the smart space to offer better and more scalable health apps. It then discussed how we envision smarter space eco-systems with a focus on one challenge: providing health app developers with tools to access local environment data. Solving this first challenge would make health apps able to perceive their environment and dynamically adjust their behavior (e.g. level of automation, need for additional information/data) to the actual status of their environment. This, in turn, would create smarter spaces and enable more scalable and reliable health apps.

In this paper, we do not address explicitly the privacy issues. However, our approach is inspired by the fog computing principles and collecting pertinent information and tagging the environment close to the users lead us to consider several questions in relation with data sovereignty/residency: What is the nature of the stored data? Where is it stored? Who holds it? Who can access it? We may consider the models [START_REF] Antignac | Privacy by Design: From Technologies to Architectures[END_REF] in the design of our tools. In our study case, we partially solve the issues by using citizen participation: Sam, the asthma suffer has the option to choose where his data comes from, where he stores it and who has access to it.

This paper envisioned smarter spaces with massive but scalable interactions between individuals and their everyday working and living environments. The potential of the smarter space is huge and not limited to health app. Smart space could help support better housing facilities, better accessibility for disabled people, smart plant, automated cooperative vehicle or car-Infrastructure cooperation.

Fig. 1 .

 1 Fig. 1. Illustration of an asthma management scenario

Fig. 3 .

 3 Fig. 3. Design of a pervasive app using cloud services/local data and local processing

Fig. 5 .

 5 Fig. 5. Building high level abstraction with raw data collected in smart environment

State 1: 78 % Raw data app Fusion tools for app developer Exposition of fused data

	Decision
	Fusion	State 1: 2 % State 2: 78 % Ignorance : 20 %
	Model	State 1: 10 % Ignorance: 90 %	State 2: 80 % Ignorance: 20 %
			Sensors


Zoïs Moitier
email: zois.moitier@univ-rennes1.fr

Stéphane Balac
email: stephane.balac@univ-rennes1.fr

Eric Darrigrand
email: eric.darrigrand-lacarrieu@univ-rennes1.fr

Monique Dauge
email: monique.dauge@univ-rennes1.fr

Yvon Lafranche
email: yvon.lafranche@univ-rennes1.fr

Fabrice Mahé
email: fabrice.mahe@univ-rennes1.fr

Moitier Zoïs

Balac Stéphane

Darrigrand Éric

Dauge Monique

Lafranche Yvon

Mahé Fabrice

Mathematical investigations of resonances in optical micro-resonators

Keywords: Scattering resonances, Optical micro-resonators, Finite Element Method, PML

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

I. INTRODUCTION

An optical micro-resonator is a device basically made of a dielectric cavity (with typical size ranging from one to several hundreds µm) coupled to waveguides or fibers for light input and output [START_REF] Heebner | Optical Microresonators: Theory, Fabrication, and Applications[END_REF]. The cavities have a characteristic set of frequencies at which it is possible to confine light-waves when the resonance conditions are met. Light passes through the cavity from the input waveguide and it builds up in intensity over multiple round-trips due to constructive interference. Depending on the optical characteristics of the cavity, it is possible to access a wide range of optical phenomena. Optical micro-resonators came to be important components in the photonic toolbox where they are implemented as lasers, amplifiers, sensors, filters, demultiplexers, switches, routers, logic gates, etc. [START_REF] Chiasera | Spherical whispering-gallery-mode microresonators[END_REF]. However, the theoretical study of optical micro-resonators in order to improve their characteristics, their design or to broaden the field of applications is currently hampered by the weakness of the mathematical models used

The authors benefit from the support of the French government "Investissements d'Avenir" program ANR-11-LABX-0020-01. This work was supported by the Interdisciplinary Mission of the French National Center for Scientific Research (CNRS) through the InFIniTi project SONATINE.

to describe their functioning and the relative inefficiency of the numerical simulation tools available.

In this work we investigate from a mathematical point of view the scattering problem associated with the computation of optical resonances in a micro-cavity for a 2D setting. The resonances are numerically computed by solving the scattering problem (detailed in Section II) by the Finite Element Method (FEM). The outgoing wave condition is handled by the use of a proper Perfectly Matched Layer (PML) presented in Section III. The properties of the PML for our scattering problem are investigated in Section IV in the particular case of a ringshaped cavity where analytical expressions for the resonances are known. We conclude with numerical simulation results in Section V.

II. THE SCATTERING PROBLEM

We consider the resonance cavity alone with invariance in the z direction. So we can assume that the optical index and the electro-magnetic field (E, H) depend only on (x, y). We denote by Ω the open connected set in R 2 filled by the dielectric cavity. We denote by ε c the dielectric permittivity of the cavity and by n c its optical index such that ε c = ε 0 n 2 c where ε 0 is the dielectric permittivity of vacuum. We also denote by Σ the boundary of Ω and by Ω the exterior domain (Ω = R 2 \ Ω) assumed to be vacuum. The optical index is the piecewise constant function n such that n = n c in Ω and n = 1 in Ω .

Under the z invariance assumption, it is well known that the set of harmonic Maxwell's equations can be split into two sub-systems of equations referred to as transverse electric (TE) equations corresponding to H = (0, 0, H z ) and transverse magnetic (TM) equations corresponding to E = (0, 0, E z ). The two subsystems being very similar, we will only consider here the (TE) case. Taking as main unknown the component u = H z of the magnetic field, the formulation of the resonance problem for TE modes reads: find

(k 0 , u) ∈ C × H 2 loc (R 2 \ Σ) with u = 0 such that div n -2 ∇u + k 2 0 u = 0 in Ω and Ω (1a) [u] = 0 and n -2 ∂ ν u = 0 across Σ (1b)
with a radiation condition at infinity that imposes that u has, outside any disk D(0, R) which contains Ω, an expansion in terms of Hankel functions of the first kind H

m of the following form in polar coordinates (r, θ):

u(r, θ) = m∈Z c m H (1) m (k 0 r) e imθ ∀θ ∈ R/2πZ. (1c)
This radiation condition describes outgoing time-harmonic waves for the time dependence e -ik0ct where c is the speed of light in vacuum. For real positive k 0 , the radiation condition (1c) coincides with the classical Sommerfeld radiation condition for outgoing waves. However it is well known, for example [START_REF] Mclean | Strongly Elliptic Systems and Boundary Integral Equations[END_REF], that problem (1) has no solution for k 0 ∈ R * and that the resonances k 0 have a negative imaginary part (Im(k 0 ) < 0) which describes the damping of the modes u in time.

From a numerical point of view, since problem (1) is set in R 2 , an artificial boundary has to be introduced to bound the computational domain. This entails the inherent difficulty of the way of choosing the boundary condition on this artificial boundary to replace the radiation condition (1c) at infinity. We have chosen to use Perfectly Matched Layers (PML) as introduced in [START_REF] Kim | The computation of resonances in open systems using a perfectly matched layer[END_REF] for the computation of resonances in open systems.

III. SCATTERING PROBLEM WITH PML

The PML method consists of a complex scaling [START_REF] Dyatlov | Mathematical theory of scattering resonances[END_REF] and a cut-off at infinity. We use a cylindrical PML with parameters σ 0 > 0, r 1 > 0 as inner radius such that Ω ⊂ D(0, r 1 ), r 2 > r 1 as outer radius, and a regularity index p ∈ {0, 1, 2, . . . }. To define the complex scaling, we introduce the function σ ∈

C p (R + ) such that σ(r) = 0 for r ≤ r 1 , σ(r) = σ 0 for r ≥ r 2 , σ(r) = σ 0 r r1 (s -r 1 ) p (r 2 -s) p ds r2 r1 (s -r 1 ) p (r 2 -s) p ds for r 1 < r < r 2 , (2)
and the following quantities:

σ(r) = σ(r) + rσ (r), d(r) = 1 + i σ(r), d(r) = 1 + i σ(r), d 0 = 1 + iσ 0 .
We set u(r, θ) = u( d(r)r, θ). The function u satisfies the system (1a,1b) in D(0, r 1 ). In the PML area for Im(d 0 k 0 ) > 0, u satisfies

∆ u + k 2 0 u = 0 for r ∈ (r 1 , +∞) \ {r 2 } (3a) [ u] = 0 and 1 d ∂ r u = 0 for r ∈ {r 1 , r 2 } (3b) lim r→+∞ u(r, θ) = 0 (3c)
where, in polar coordinates (r, θ),

∆ u = 1 r d ∂ ∂r r d d ∂ u ∂r + 1 d 2 r 2 ∂ 2 u ∂θ 2 .
Choosing r 3 ≥ r 2 , we approach u by u that satisfies the system (1a,1b) in D(0, r 1 ) and

∆ u + k 2 0 u = 0 for r ∈ (r 1 , r 3 ) \ {r 2 } (4a) [ u] = 0 and 1 d ∂ r u = 0 for r ∈ {r 1 , r 2 } (4b) u = 0 for r = r 3 (4c)
The weak formulation of the eigenvalue problem set in the bounded domain D(0, r 3 ) reads: find ( k 0 , u) ∈ C × H 1 0 (D(0, r 3 )) with u = 0 and Im( k 0 ) < 0 such that ∀v ∈ H 1 0 (D(0, r 3 ))

D(0,r3) n -2 M ∇ u • ∇v dxdy -k 2 0 D(0,r3) d d u v dxdy = 0 (5)
where M is the matrix

M = 1 x 2 + y 2   x 2 d d + y 2 d d xy d d -d d xy d d -d d x 2 d d + y 2 d d  
Problem (5) can then be solved by the Finite Element Method in the usual way. For the present study, we have used the open-source Finite Element Library XLIFE++ [START_REF]XLiFE++ : an extended library of finite element library in C++[END_REF] IV. A REFERENCE PROBLEM

In order to investigate computational issues when solving problem (1) by the FEM using PML, we consider a "reference problem" where the micro-cavity has a ring shape with inner radius r i and outer radius r o . For such a geometry, separating the variables in polar coordinates provides analytical expressions for the optical modes. We find that solutions u to problem (1) can be expressed as u(r, θ) = m∈Z u m (r) e imθ where u m is given in terms of Bessel's functions J m , Y m , and

H (1) m as u m (r) =    C 1 J m (k 0 r) if r < r i C 2 J m (n c k 0 r) + C 3 Y m (n c k 0 r) if r ∈ (r i , r o ) C 4 H (1) m (k 0 r) if r > r o
where the four constants C 1 , . . . , C 4 are determined by the interface conditions (1b) and the resonance k 0 is obtained by solving the "modal equation"

J m (nck0ri) nc J m (k0ri) -Jm(nck0ri) Jm(k0ri) Y m (nck0ri) nc J m (k0ri) -Ym(nck0ri) Jm(k0ri) = Jm(nck0ro) H (1) m (k0ro) - J m (nck0ro) nc H (1) m (k0ro) Ym(nck0ro) H (1) m (k0ro) - Y m (nck0ro) nc H (1) m (k0ro) (6)
If m∈Z u m (r) e imθ is a resonance for k 0 then each term u m (r) e imθ is a resonance for the same k 0 . We have conducted numerical experiments on this reference problem to explore the influence of the different parameters involved in the definition of the PML on the accuracy of the computations.

First, we have noted that a much better approximation is obtained by increasing the FE polynomial order rather than increasing the number of elements. We have also observed that the parameter σ 0 in (2) must be chosen carefully: large enough to ensure a good convergence of the approximation, but not too large in which case the convergence rate would be affected. Finally, a better performance is obtained by using the lowest regularity on σ corresponding to p = 0 in (2). However, choosing p = 0 and σ ∈ C 0 ([0, r 2 ]) requires that the mesh of the domain respects the discontinuity of σ at r = r 1 . We can choose r 3 = r 2 , this is simpler to code and does not alter the convergence. Other experiments show that the PML width r 2 -r 1 does not need to be large.

V. NUMERICAL SIMULATION RESULTS

We present numerical simulations obtained for a ring with inner radius r i = 0.3, outer radius r o = 0.5, and optical index n c = 10. Similar simulations (on disk) can be found in [START_REF] Nannen | Spurious modes of the complex scaled Helmholtz equation[END_REF], [START_REF] Moiola | Acoustic transmission problems: wavenumber-explicit bounds and resonance-free regions[END_REF]. We compare a 1D calculation (done with FE of order 30 and around 200 dofs-as precise as the analytic formulas presented in the previous paragraph) and a 2D calculation (done with triangular FE of order 6 and around 60000 dofs) of the resonances. Computations were carried out using XLIFE++ [START_REF]XLiFE++ : an extended library of finite element library in C++[END_REF] with the following PML parameters: r 1 = 1.1, r 2 = r 3 = 1.6, σ 0 = 5, p = 0. We show on Fig. 1 the optical modes corresponding to the following resonances k 0 : j m 1D 2D rel. diff. 1 5 2.000254 -i 10 -7.859 2.000254 -i 10 -7.726 10 -6.5 1 10 2.942771 -i 10 -14.09 2.942771 -i 10 -10.64 10 -7.4 1 15 4.004424 -i 10 -14.79 4.004424 -i 10 -13.79 10 -7.3 2 5 3.367404 -i 10 -5.329 3.367406 -i 10 -5.336 10 -6.2 2 10 4.050435 -i 10 -13.44 4.050437 -i 10 -10.23 10 -6.3 2 15 4.975665 -i 10 -16.82 4.975673 -i 10 -14.14 10 -5.7 3 5 4.842406 -i 10 -3.657 4.842411 -i 10 -3.670 10 -5.7 3 10 5.354871 -i 10 -10.89 5.354881 -i 10 -9.789 10 -5.7 3 15 6.095389 -i 10 -14.27 6.095450 -i 10 -11.60 10 -5.0 We have also considered a race-track shaped cavity with the following features: length of the straight line 1, external turn radius 0.5, width of the cavity 0.2, and optical index n c = 10.

VI. CONCLUSION

We have set up a way to compute resonances in 2D dielectric micro-cavity, so we can find resonances on geometries where no analytical solution is available. Furthermore, we know [START_REF] Lam | Explicit asymptotic formulas for the positions, widths, and strengths of resonances in Mie scattering[END_REF] that the resonances have asymptotic expansions as m tends to infinity for a disk. For our computations on the ring, the resonances satisfy the same asymptotic. In a future work, we plan to investigate if there are asymptotic expansions for other geometries.

Fig. 1 .

 1 Fig. 1. Optical modes for a ring cavity (real part of the eigenfunction u) for various resonances k m 0,j corresponding to j = 1, 2, 3 (row index) and m = 5, 10, 15 (column index).

  Fig 2 shows the optical mode corresponding to k 0 = 3.19-i 10 -4.13 . Computations were carried out with triangular FE of order 6, around 60000 dofs, and the following PML parameters: r 1 = 1.1, r 2 = r 3 = 1.6, σ 0 = 5, p = 0.

Fig. 2 .

 2 Fig. 2. Optical modes for a race-track shaped cavity (real part of the eigenfunction u) for k 0 = 3.19 -i 10 -4.13 .

ACKNOWLEDGMENT

The authors thank Y. Dumeige, P. Féron and Y. Leger from FOTON Laboratory (CNRS UMR 6082) for their contributions to the SONATINE project.