F Dubois
email: dubois@lmgc.univ-montp2.fr

M Jean

Frictional Contact Numerical Models for Numerous Collections of Rigid or Deformable Bodies

Interactions between contacting bodies may be physically quite complicated. Deriving a revelant numerical model for dynamical problem raises many issues especially due to the fact that it mixes various time and space scales. The numerical modelling, i.e. the writing of an interaction law and the choice of a numerical strategy, depends on the problem into consideration and what are the main results expected from a numerical simulation. Primary phenomena are impenetrability of materials and dry friction. Other phenomena like cohesive behaviour, local deformability, plastic behaviour, viscous friction, wear, etc. should be taken into account. Large collections of bodies, like granular materials, have their own specificities :

• the number of contacts may be very large, 1000 for small samples to 300.000 for "representative" samples;

• dynamics is to be taken into account since the motion equations governing the grains suppose they are free to move; meanwhile they are submitted to stiff interaction laws;

• simultaneous contacts are occuring;

• the solutions may suffer from indetermination; in rigid bodies collections, an equilibrium state might be realized by an infinity of combinations of reactions belonging to some subset. Coulomb's law is responsible for this non uniqueness.

Roughly, there are two numerical strategies :

• (SCD) The interaction laws are described by mappings relating the main local unknowns involved in the contact -the gap, the relative velocity, and the dual variable, the reaction (a force or an impulse). One ends up with a stiff ordinary differential equation which may be treated with an explicit integration scheme. Small time steps are necessary, together with some kind of damping ensuring the scheme stability. The episods of contact are finely described. The so called Distinct Element Method (P. Cundall & al.) or the Molecular Dynamics methods belong to this Smooth Contact Dynamics (SCD) class.

• (NSCD) The interaction laws are described by Signorini condition and Coulomb's law, or other so called "derived laws". These laws are written as multi mappings relating the local unknowns, in the frame of Convex Analysis following J.J. Moreau, together with first order schemes for the dynamical equations. In the NSCD methods, at each contact, the impulse over the time step and the relative velocity at the end of the time step are computed together with the dynamical equation terms so as to satisfy Signorini condition and Coulomb's law. These impulses account for all events supposed to occur during a time step, free flights, action of external forces, single or numerous impacts between pairs or agglomerates of contacting bodies. Iterations must be performed in order to exhibit the values of the unknowns. The method allows large time steps. This is the Non Smooth Contact Dynamics (NSCD), a fully implicit method, initiated by J.J. Moreau [START_REF] Moreau | Some basics of unilateral dynamics[END_REF], M. Jean [START_REF] Jean | Engineering, special issue on computational modeling of contact and friction[END_REF] and implemented in the LMGC90 software by F. Dubois and M. Jean.

Several elementary examples illustrating the behaviour of the SCD and the NSCD method will be given concerning the impact of a few particles, using non smooth laws or spring damper contacting devices, using SCD or NSCD methods with small or large time steps: a bouncing ball, billiard balls, a chain of cohesive balls, etc. It appears that indeed, explicit SCD methods give a fine description of the contact episod to be computed with small time steps. As far as large collections of bodies are concerned, coming into details of episods of contact might be irrelevant or superfluous. For instance SCD methods generates genuine or artefact propagation of waves which are not of interest in a soil mechanic quasistatic experiment. Furthermore the mechanical values of the spring damper devices are hardly known. These devices should be understood as a smoothing or penalization process, the harder the springs, the smaller the time steps. The implicit NSCD method uses large time steps and iterations are necessary in order to select a solution. The question of indetermination i.e. the fact that an equilibrium state of a sample may be realized by an infinity of combinations of reactions belonging to some subset, is raised in full rigid bodies models. This is illustrated by the example of a wedged disk and of the billiard balls. When the Signorini Coulomb interaction law is supplemented with a spring damper device, in the spirit of derived Signorini Coulomb's law (a smoother alternative way to comply with unilaterality) the indetermination is partly suppressed. Unfortunatly it appears that those models suffer from sensitivity with respect to smoothing parameters, which exhibits a behaviour similar to indetermination. It seems that as far one is concerned with average data as stresses, resulting pressures on walls, or statistic distributions of forces, these indeterminations do not matter much because they hit mainly the weak forces.


Timothée Gerber

Martin Dutasta

Laurent Girin

Cédric Févotte
email: cedric.fevotte@telecom-paristech.fr

PROFESSIONALLY-PRODUCED MUSIC SEPARATION GUIDED BY COVERS

This paper addresses the problem of demixing professionally produced music, i.e., recovering the musical source signals that compose a (2-channel stereo) commercial mix signal. Inspired by previous studies using MIDI synthesized or hummed signals as external references, we propose to use the multitrack signals of a cover interpretation to guide the separation process with a relevant initialization. This process is carried out within the framework of the multichannel convolutive NMF model and associated EM/MU estimation algorithms. Although subject to the limitations of the convolutive assumption, our experiments confirm the potential of using multitrack cover signals for source separation of commercial music.

INTRODUCTION

In this paper, we address the problem of source separation within the framework of professionally-produced (2channel stereo) music signals. This task consists of recovering the individual signals produced by the different instruments and voices that compose the mix signal. This would offer new perspectives for music active listening, editing and post-production from usual stereo formats (e.g., 5.1 upmixing), whereas those features are currently roughly limited to multitrack formats, in which a very limited number of original commercial songs are distributed.

Demixing professionally produced music (PPM) is particularly difficult for several reasons [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF][START_REF] Ozerov | Multichannel nonnegative tensor factorization with structured constraints for user-guided audio source separation[END_REF][START_REF] Sturmel | Linear mixing models for active listening of music productions in realistic studio condition[END_REF]. Firstly, the mix signals are generally underdetermined, i.e., there are more sources than mix channels. Secondly, some sources do not follow the point source assumption that is often implicit in the (convolutive) source separation models of the signal processing literature. Also, some sources can be panned in the same direction, convolved with large reverberation, or processed with artificial audio effects that are more or less easy to take into account in a separation framework. PPM separation is thus an ill-posed problem and separation methods have evolved from blind to informed source separation (ISS), i.e., methods that exploit Permission to make digital or hard copies of all or part of this work for personal or classroom use is granted without fee provided that copies are not made or distributed for profit or commercial advantage and that copies bear this notice and the full citation on the first page. c 2012 International Society for Music Information Retrieval. some "grounded" additional information on the source/mix signals and mix process. For example, the methods in [START_REF] Dubnov | Optimal filtering of an instrument sound in a mixed recording using harmonic model and score alignment[END_REF][START_REF] Ewert | Score-informed voice separation for piano recordings[END_REF][START_REF] Ewert | Using score-informed constraints for NMF-based source separation[END_REF][START_REF] Ganseman | Source separation by score synthesis[END_REF][START_REF] Woodruff | Remixing stereo music with score-informed source separation[END_REF] exploit the musical score of the instrument to extract sources, either directly or through MIDI signal synthesis. In user-guided approaches, the listener can assist the separation process in different ways, e.g., by humming the source to be extracted [START_REF] Smaragdis | Separation by "humming": User-guided sound extraction from monophonic mixtures[END_REF], or by providing information on the sources direction [START_REF] Vinyes | Demixing commercial music productions via human-assisted timefrequency masking[END_REF] or temporal activity [START_REF] Ozerov | Multichannel nonnegative tensor factorization with structured constraints for user-guided audio source separation[END_REF]. An extreme form of ISS can be found in [START_REF] Faller | Enhancing stereo audio with remix capability[END_REF][START_REF] Gorlow | Informed source separation: Underdetermined source signal recovery from an instantaneous stereo mixture[END_REF][START_REF] Liutkus | Informed source separation through spectrogram coding and data embedding[END_REF][START_REF] Parvaix | Informed source separation of linear instantaneous under-determined audio mixtures by source index embedding[END_REF][START_REF] Parvaix | A watermarking-based method for informed source separation of audio signals with a single sensor[END_REF] and in the Spatial Audio Object Coding (SAOC) technology recently standardized by MPEG [START_REF] Engdegård | MPEG spatial audio object coding-the ISO/MPEG standard for efficient coding of interactive audio scenes[END_REF]: here, the source signals themselves are used for separation, which makes sense only in a coder-decoder configuration.

In the present paper, we remain in the usual configuration where the original multitrack signals are not available, although we keep the latter spirit of using source signals to help the demixing process: we propose to use cover multitrack signals for this task. This idea is settled on several facts. Firstly, a cover song can be quite different from the original for the sake of artistic challenge. But very interestingly, for some applications/markets a cover song is on the contrary intended to be as close as possible to the original song: instruments composition and color, song structure (chorus, verses, solos), and artists interpretation (including the voices) are then closely fitted to the original source signals, hence having a potential for source separation of original mixes. Remarkably, it happens that multitracks of such "mimic" covers are relatively easy to find on the market for a large set of famous pop songs. In fact, they are much easier to obtain than original multitracks. This is because the music industry is very reluctant to release original works while it authorizes the licensed production of mimic multitracks on a large scale. In the present study, we will use such multitracks provided by iKlax Media which is a partner of the DReaM project. 1 iKlax Media produces software solutions for music active listening and has licensed the exploitation of a very large set of cover multitracks of popular songs. Therefore, this work involves a sizeable artistic and commercial stake. Note that similar material can be obtained from several other companies.

We set the cover-informed source separation principle within the currently very popular framework of separation methods based on a local time-frequency (TF) complex Gaussian model combined with a non-negative matrix factorization (NMF) model for the source variances [START_REF] Févotte | Nonnegative matrix factorization with the Itakura-Saito divergence. With application to music analysis[END_REF][START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF][START_REF] Ozerov | A general flexible framework for the handling of prior information in audio source separation[END_REF].

Iterative NMF algorithms for source modeling and separation have shown to be very sensitive to initialization. We turn this weakness into strength within the following twostep process in the same spirit as the work carried out on signals synthesized from MIDI scores in, e.g., [START_REF] Ganseman | Source separation by score synthesis[END_REF] or by humming in [START_REF] Smaragdis | Separation by "humming": User-guided sound extraction from monophonic mixtures[END_REF]. First, source-wise NMF modeling is applied on the cover multitrack, and the result is assumed to be a suitable initialization of the NMF parameters of the original sources (that were used to produce the commercial mix signal). Starting from those initial values, the NMF process is then refined by applying to the mix the convolutive multichannel NMF model of [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF]. This latter model provides both refined estimation of the source-within-mix (aka source images) NMF parameters and source separation using Wiener filters built from those parameters.

The paper is organized as follows. In Sections 2 and 3, we respectively present the models and method employed. In Sections 4 and 5, we present the experiments we conducted to assess the proposed method, and in Section 6, we address some general perspectives.

FRAMEWORK: THE CONVOLUTIVE MULTICHANNEL NMF MODEL

Mixing Model

Following the framework of [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF], the PPM multichannel mix signal x(t) is modeled as a convolutive noisy mixture of J source signals s j (t). Using the short-time Fourier transform (STFT), the mix signal is approximated in the TF domain as:

x f n = A f s f n + b f n , (1)
where x f n = [x 1,f n , . . . , x I,f n ] T is the vector of complexvalued STFT coefficients of the mix signal, s f n = [s 1,f n , . . . , s J,f n ] T is the vector of complex-valued STFT coefficients of the sources,

b f n = [b 1,f n , . . . , b I,f n ] T is a zero- mean Gaussian residual noise, A f = [a 1,f , . . . , a J,f
] is the frequency-dependent mixing matrix of size I × J (a j,f is the mixing vector for source j), f ∈ [0, F -1] is the frequency bin index and n ∈ [0, N -1] is the time frame index. This approach implies standard narrowband assumption (i.e., the time-domain mixing filters are shorter than the STFT window size).

Source model

Each source s j,f n is modeled as the sum of

K j latent com- ponents c k,f n , k ∈ K j , i.e., s j,f n = k∈Kj c k,f n , (2)
where {K j } j is a non-trivial partition of {1, . . . , K}, K ≥ J (K j is thus the cardinal of K j ). Each component c k,f n is assumed to follow a zero-mean proper complex Gaussian distribution of variance w f k h kn , where

w f k , h kn ∈ R + , i.e., c k,f n ∼ N c (0, w f k h kn ).
The components are assumed to be mutually independent and individually inde-pendent across frequency and time, so that we have:

s j,f n ∼ N c (0, k∈Kj w f k h kn ). (3)
This source model corresponds to the popular non-negative matrix factorization (NMF) model as applied to the source power spectrogram

|S j | 2 = {|s j,f n | 2 } f n : |S j | 2 ≃ W j H j , (4)
with non-negative matrices W j = {w f k } f,k∈Kj of size F × K j and H j = {h kn } k∈Kj ,n of size K j × N . The columns of W j are generally referred to as spectral pattern vectors, and the rows of H j are referred to as temporal activation vectors. NMF is largely used in audio source separation since it appropriately models a large range of musical sounds by providing harmonic patterns as well as non-harmonic ones (e.g., subband noise).

Parameter estimation and source separation

In the source modeling context, the NMF parameters of a given source signal can be obtained from the observation of its power spectrogram using Expectation-Maximization (EM) iterative algorithms [START_REF] Févotte | Nonnegative matrix factorization with the Itakura-Saito divergence. With application to music analysis[END_REF]. In [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF], this has been generalized to the joint estimation of the J sets of NMF source parameters and I × J × F mixing filters parameters from the observation of the mix signal power spectrogram. More precisely, two algorithms were proposed in [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF]. An EM algorithm consists of maximizing the exact joint likelihood of the multichannel data, whereas a multiplicative updates (MU) algorithm, maximizes the sum of individual channel log-likelihood. If the former better exploits the interchannel dependencies and gives better separation results,2 the latter has a lower computation cost. Those algorithms will not be described in the present paper, the reader is referred to [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF] for technical details. Once all the parameters are estimated, the source signals (or their spatial images y j,f n = a j,f s j,f n ) are estimated using spatial Wiener filtering of the mix signal:

ŝfn = Σ s,f n A H f Σ -1 x,f n x f n , (5)
where Σ s,f n is the (estimated) covariance matrix of the source signals, and

Σ x,f n = A f Σ s,f n A H f + Σ b,f
is the (estimated) covariance matrix of the mix signal.

PROPOSED COVER-INFORMED SEPARATION TECHNIQUE

Cover-based initialization

It is well-known that NMF decomposition algorithms are highly dependent on the initialization. In fact, the NMF model does not guarantee the convergence to a global minimum but only to a local minimum of the cost function, making a suitable initialization crucial for the separation performance. In the present study, we have at our disposal the 2-channel stereo multitrack cover of each song to separate, and the basic principle is to use the cover source tracks to provide relevant initialization for the joint multichannel decomposition. Therefore, the NMF algorithms mentioned in Section 2 are applied on PPM within the following configuration. A first multichannel NMF decomposition is run on each stereo source of the cover multitrack (with random initialization). Thus, we obtain a modeled version of each cover source signal in the form of three matrices per source: W cover j , H cover j and A cover j = {a cover ij,f } i∈ [START_REF] Dubnov | Optimal filtering of an instrument sound in a mixed recording using harmonic model and score alignment[END_REF][START_REF] Duong | Underdetermined reverberant audio source separation using a full-rank spatial covariance model[END_REF],f . The results are ordered according to:

W mix init = [W cover 1 . . . W cover J ] (6)
H mix init =    H cover 1 . . . H cover J    (7)
A mix init = [A cover 1 . . . A cover J ] (8)
Then, ( 6), [START_REF] Févotte | Nonnegative matrix factorization with the Itakura-Saito divergence. With application to music analysis[END_REF], and ( 8) are used as an initialization for a second convolutive stereo NMF decomposition run on the mix signal as in [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF]. During this second phase, the spectral pattern vectors and time activation vectors learned from the cover source tracks are expected to evolve to match the ones corresponding to the signals used to produce the commercial mix, while the resulting mixing vectors are expected to fairly model the mix process.

Pre-processing: time alignment of the cover tracks

One main difference between two versions of the same music piece is often the temporal misalignment due to both tempo variation (global misalignment) and musical interpretation (local misalignments). In a general manner, time misalignment can corrupt the separation performances if the spectral pattern vectors used for initialization are not aligned with the spectral patterns of the sources within the mix. In the present framework, this problem is expected to be limited by the intrinsic automatic matching of temporal activity vectors within the multichannel NMF decomposition algorithm. However, the better the initial alignment, the better the initialization process and thus expected final result. Therefore, we limit this problem by resynchronizing the cover tracks with the mix signal, in the same spirit as the MIDI score-to-audio alignment of [START_REF] Ewert | Using score-informed constraints for NMF-based source separation[END_REF] or the Dynamic Time Warping (DTW) applied on synthesized signals in [START_REF] Ganseman | Source separation by score synthesis[END_REF]. In the present study, this task is performed at quarter-note accuracy using the Beat Detective tool from the professional audio editing software Avid ProTools R . This step allows minimizing synchronization error down to less than a few TF frames, which is in most cases below the synchronization error limit of 200 ms observed in [START_REF] Ewert | Using score-informed constraints for NMF-based source separation[END_REF].

In-depth study of desynchronization on source separation is kept for future works.

Exploiting the temporal structure of source signals

In order to further improve the results, we follow a userguided approach as in [START_REF] Ozerov | Multichannel nonnegative tensor factorization with structured constraints for user-guided audio source separation[END_REF]. The coefficients of matrix H are zeroed when the source is not active in the mix, exploiting audio markers of silence zones in the cover source tracks. As there still may be some residual misalignment between the commercial song and the cover after the preprocessing, we relax these constraints to 3 frames before and after the active zone. When using the MU algorithm, the zeroed coefficients remain at zero. When using the EM algorithm, the update rules do not allow the coefficients of H to be strictly null, hence, we set these coefficients to the eps value in our Matlab R implementation. Observations confirm that these coefficients remain small throughout all the decomposition.

Summarizing the novelty of the proposed study

While our process is similar in spirit to several existing studies, e.g., [START_REF] Ewert | Using score-informed constraints for NMF-based source separation[END_REF][START_REF] Ganseman | Source separation by score synthesis[END_REF][START_REF] Smaragdis | Separation by "humming": User-guided sound extraction from monophonic mixtures[END_REF], our contribution to the field involves:

• the use of cover multitrack signals instead of hummed or MIDI-synthesis source signals. Our cover signals are expected to provide a more faithful image of the original source signals in the PPM context.

• a stereo NMF framework instead of a mono one. The multichannel framework is expected to exploit spatial information in the demixing process (as far as the convolutive model is a fair approximation of the mixing process). It provides optimal spatial Wiener filters for the separation, as opposed to the {estimated magnitude + mix phase} resynthesis of [START_REF] Ganseman | Source separation by score synthesis[END_REF] or the (monochannel) soft masks of [START_REF] Smaragdis | Separation by "humming": User-guided sound extraction from monophonic mixtures[END_REF].

• a synchronization pre-process relying on tempo and musical interpretation instead of, e.g., frame-wise DTW. This is completed with the exploitation of the sources temporal activity for the initialization of H.

EXPERIMENTS

Data and experimental settings

Assessing the performances of source separation on true professionally-produced music data is challenging since the original multitrack signals are necessary to perform objective evaluation but they are seldom available. Therefore, we considered the following data and methodology. The proposed separation algorithm was applied on a series of 4 well-known pop-music songs for which we have the stereo commercial mix signal and two different stereo multitrack covers (see Table 2). The first multitrack cover C1 was provided by iKlax Media, and the second one C2 has been downloaded from the commercial website of another company. We present two testing configurations:

• Setting 1: This setting is used to derive objective measures (see below). C1 is considered as the "original multitrack", and used to make a stereo remix of the song which is used as the target mix to be separated. This remix has been processed by a qualified sound engineer with a 10-year background in music • Setting 2: The original commercial mix is separated using C1 as the cover. This setting is used for subjective evaluation in real-world configuration.

The covers are usually composed of 8 tracks which are quite faithful to the commercial song content as explained in the introduction. For simplicity we merged the tracks to obtain 4 to 6 source signals. 4 All signals are resampled at 32kHz, since source separation above 16kHz has very poor influence on the quality of separated signals and this enables to reduce computations. The experiments are carried out on 30s excerpts of each song.

It is difficult to evaluate the proposed method in reference to existing source separation methods since the cover information is very specific. However, in order to have a reference, we also applied the algorithm with a partial initialization: the spectral patterns W are here initialized with the cover spectral patterns, whereas the time activation vectors H are randomly initialized (vs. NMF initialization in the full cover-informed configuration). This enables to i) separate the contribution of cover temporal information, and ii) simulate a configuration where a dictionary of spectral bases is provided by an external database of instruments and voices. This was performed for both EM and MU algorithms. The main technical experimental parameters are summarized in Table 1.

Separation measures

To assess the separation performances in Setting 1, we computed the signal-to-distortion ratio (SDR), signal-tointerference ratio (SIR), signal-to-artifact ratio (SAR) and source image-to-spatial distortion ratio (ISR) defined in [START_REF] Vincent | Performance measurement in blind audio source separation[END_REF]. We also calculated the input SIR (SIR in ) defined as the ratio between the power of the considered source and 3 The source images are here the processed version of C1 just before final summation, hence we do not consider post-summation (non-linear) processing. The consideration of such processing in ISS, as in, e.g., [START_REF] Sturmel | Linear mixing models for active listening of music productions in realistic studio condition[END_REF], is part of our current efforts. 4 The gathering was made according to coherent musical sense and panning, e.g., grouping two electric guitars with the same panning in a single track. It is necessary to have the same number of tracks between an original version and its cover. Furthermore, original and cover sources should share approximately the same spatial position (e.g., a cover version of a left panned instrument should not be right panned!) the power of all the other sources in the mix to be separated. We consider this criterion because all sources do not contribute to the mix with the same power. Hence, a source with high SIR in is easier to extract than a source with a low SIR in , and SIR in is used to characterize this difficulty.

RESULTS

Objective evaluation

Let us first consider the results obtained with Setting 1. The results averaged across all sources and songs are provided in Table 3. The maximal average separation performance is obtained with the EM cover-informed algorithm with SDR = 2.45dB and SIR = 4.00dB. This corresponds to a source enhancement of SDR -SIR in = 10.05dB and SIR-SIR in = 11.60dB, with the average global SIR in being equal to -7.60dB. These results show that the overall process leads to fairly good source reconstruction and rejection of competing sources. Figure 1a illustrates the separation performances in terms of the difference SDR -SIR in for the song "I will survive". The separation is very satisfying for tracks with sparse temporal activity such as Brass. The Strings track, for which the point source assumption is less relevant, obtains correct results, but tends to spread over other sources images such as Bass. Finally, when cover tracks musically differ from their original sources, the separation performance decreases. This is illustrated with the Electric Guitar (EGtr) and Bass tracks, which do not fully match the original interpretation.

Let us now discuss the cover informed EM and MU methods in relation to the initialization of spectral bases only, referred to as W init . The cover-based EM algorithm provides a notable average SDR improvement of 2.41dB over EM with W init initialization, and a quite large improvement in terms of SIR (+5.97dB), hence a much better interference rejection. The cover-based MU algorithm also outperforms the MU W init configuration to the same extent (e.g., +2.36dB SDR and +6.18dB SIR improvement). This reveals the ability of the method to exploit not only spectral but also temporal information provided by covers.

Note that both cover-based and W init EM methods outperform the corresponding MU methods in terms of SDR. However, it is difficult to claim for clear-cut EM's better use of the inter-channel mutual information, since EM is slightly lower than MU for SIR (approx. 4dB vs. 5dB for the cover-informed method). In fact, the multichannel framework can take advantage of both spectral and spatial information for source extraction, but this depends on the source properties and mixing configuration. In the song "Walk this way", which detailed results are given in Figure 1b, all sources but the Electric Guitar 1 (Egtr1) are panned at the center of the stereo mixture. Thus, the SDR -SIR in obtained for Egtr1 reaches 20.32dB, as the algorithm relies strongly on spatial information to improve the separation. On the other hand, the estimated Vocal track in "I will survive" is well separated (+8.57dB SDR -SIR in for the cover-informed EM) despite being centered and coincident to other tracks such as Bass, Drums and Electric Guitar (EGtr). In this case, the proposed multichannel NMF framework seems to allow separation of spatially coincident sources with distinct spectral patterns. Depending on the song, some sources obtain better SDR results with the MU algorithm. For example, in "Walk this way", the SDR -SIR in for the Drums track increased from 6.59dB with the EM method to 9.74dB with the MU method. As pointed out in [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF], the point source assumption certainly does not hold in this case. The different elements of the drums are distributed between both stereo channels and the source image cannot be modeled efficiently as a convolution of a single point source. By discarding a large part of the inter-channel information, the MU algorithm gives better results in this case. Preliminary tests using a monochannel NMF version of the entire algorithm (monochannel separation using monochannel initialization, as in, e.g., [START_REF] Ganseman | Source separation by score synthesis[END_REF][START_REF] Smaragdis | Separation by "humming": User-guided sound extraction from monophonic mixtures[END_REF]), even show slightly better results for the Drums track, confirming the irrelevancy of the point source convolutive model in this case.

Finally, it can be mentioned that the number of NMF components per source K j does not influence significantly the SDR and SIR values, although we perceive a slight improvement during subjective evaluation for K j = 50. 5

Discussion

Informal listening tests on the excerpts from Setting 2 confirm the previous results and show the potential of coverinformed methods for commercial mix signal separation. 6 Our method gives encouraging results on PPM when point 5 Assessing the optimal number of components for each source is a challenging problem left for future work. 6 Examples of original and separated signals are available at http://www.gipsa-lab.grenoble-inp.fr/∼laurent.girin/demo/ismir2012.html. source and convolutive assumptions are respected. For instance, the vocals are in most cases suitably separated, with only long reverberation interferences. As expected, the quality of the mix separation relies on the quality and faithfulness of the cover. A good point is that when original and cover interpretations are well matched, the separated signal sounds closer to the original than to the cover, revealing the ability of the adapted Wiener filters to well preserve the original information.

Comparative experiments with spectral basis initialization only (W init ) confirm the importance of the temporal information provided by covers, Although this has not been tested formally, the cover-to-mix alignment of Section 3.2 was shown by informal tests to also contribute to good separation performances.

CONCLUSION

The results obtained by plugging the cover-informed source separation concept in the framework of [START_REF] Ozerov | Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation[END_REF] show that both spectral and temporal information provided by cover signals can be exploited for source separation. This study indicates the interest (and necessity) of using high-quality covers. In this case, the separation process may better take into consideration the music production subtleties, compared to MIDI-or hummed-informed techniques.

Part of the results show the limitations of the convolutive mixing model in the case of PPM. This is the case for sources that cannot be modeled efficiently as a point source convolved on each channel with a linear filter, such as large instruments (e.g., drums and piano). Also, some tracks such as vocals make use of reverberation times much higher than our analysis frame. As a result, most of the vocals reverberation is not properly separated. The present study and model also do not consider the possible nonlinear processes applied during the mixing process.

Therefore, further research directions include the use of more general models for both sources and spatial processing. For instance, we plan to test the full-rank spatial covariance model of [START_REF] Duong | Underdetermined reverberant audio source separation using a full-rank spatial covariance model[END_REF], within the very recently proposed general framework of [START_REF] Ozerov | A general flexible framework for the handling of prior information in audio source separation[END_REF] which also enables more specific source modeling, still in the NMF framework (e.g., source-filter models). Within such general model, sources actually composed of several instruments (e.g., drums) may be spectrally and spatially decomposed more efficiently and thus better separated.

Figure 1 :

 1 Figure 1: Separation results

Table 1 :

 1 Experimental settings production, using Avid ProTools R .3 C2 is considered as the cover version and is used to separate the target mix made with C1.

	Tracks duration	30 s
	Number of channels	I=2
	Sampling Rate	32 kHz
	STFT frame size	2048
	STFT overlap	50 %
	Number of iterations	500
	Number of NMF components 12 or 50

Table 2 :

 2 Experimental dataset

	Method	SDR ISR	SIR	SAR
	EM W init	0,04 3,51 -1,96 4,82
	EM Cover-based	2.45 6.58 4.00	5.38
	EM Improvement 2,41 3,08 5,97	0,56
	MU W init	-0,98 3,58 -1,14 3,40
	MU Cover-based	1.38 6.83 5.04	2.95
	MU Improvement 2,36 3,24 6,18 -0,45

Table 3 :

 3 Average source separation performance for 4 PPM mixtures of 4 to 6 sources (dB).

This research is partly funded by the French National Research Agency (ANR) -Grant CONTINT 09-CORD-006.

When point source and convolutive mixing assumptions are verified.


A Parret-Freaud

P Gosselet
email: gosselet@lmt.ens-cachan.fr

C Rey
email: rey@lmt.ens-cachan.fr

CONTROL OF DISCRETIZATION ERROR AND ADAPTATIVITY FOR FINITE ELEMENT PROBLEMS SOLVED BY DOMAIN DECOMPOSITION

Keywords: A posteriori error estimation, verification, adaptativity, non-overlapping domain decomposition, FETI

This paper presents a strategy for adaptative verification of computations driven in a sub-structured framework arising from the use of non-overlapping domain decomposition methods. We currently restrict ourselves to global estimates of discretization error in the context of linear mechanical problems. The key point of our method which is based on the use of constitutive relation error is a parallel recovery of admissible interface fields. This allows to obtain relevant estimates far before the domain decomposition solver reaches convergence.

INTRODUCTION

For the last decades, three trends have grown and reinforced each other in computational mechanics: the fast growth of hardware computational capacities, the requirement of finer and larger finite element models for industrial simulations and the development of efficient multiscale computational strategies amongst which non-overlapping domain decomposition (DD) methods (like FETI [START_REF] Farhat | Implicit parallel processing in structural mechanics[END_REF] and BDD [START_REF] Le | Domain decomposition methods in computational mechanics[END_REF][START_REF] Mandel | Balancing domain decomposition[END_REF]) are very popular since they have proved to be scalable in many applications. Unfortunately, the quality of the discretized models is not always checked in spite of the development of verification techniques and error estimators [START_REF] Ladeveze | Error estimate procedure in the finite element method and applications[END_REF][START_REF] Babuška | Error estimates for adaptive finite element computation[END_REF][START_REF] Pares | Subdomain-based flux-free a posteriori error estimators[END_REF][START_REF] Zienkiewicz | The superconvergent patch recovery and a posteriori error estimates. part 2: Error estimates and adaptivity[END_REF].

Our objective is to warranty the quality of computations of large finite element models by providing a strategy to control discretization error in large FE problems solved by non-overlapping domain decomposition methods.

We first focus on global a posteriori error estimates and associated element-wise contributions obtained through the error in constitutive relation [START_REF] Ladeveze | Error estimate procedure in the finite element method and applications[END_REF] for linear mechanical problems. Our method [START_REF] Parret-Fréaud | Fast estimation of discretization error for fe problems solved by domain decomposition[END_REF] has been developed within a generic framework for domain decomposition [START_REF] Gosselet | Non-overlapping domain decomposition methods in structural mechanics[END_REF]. It computes error estimate in a fully parallel way for both primal (BDD) and dual (FETI) approaches whatever the state (converged or not) of the associated iterative solver. To do so, our strategy consists in building, in parallel and during the iterations, displacement and stress fields which are kinematically admissible (KA) and statically admissible (SA) on the whole structure. The numerical costs remain low thanks to the exploitation of classical domain decomposition preconditioners to build in parallel continuous displacement and balanced traction fields independently on each subdomain. It comes out that the study convergence of the discretization error estimator with respect to the interface error of the domain decomposition method (displacement gap and lack of balance of interface forces) enables to propose new stopping criteria for the DD solver.

Though our method has been developed for both dual and primal approach of domain decomposition, this paper aims to present our strategy especially in the framework of the FETI method. We thus invite the interested reader to refer to [START_REF] Parret-Fréaud | Fast estimation of discretization error for fe problems solved by domain decomposition[END_REF] for more details about the specific setup related to the BDD method. Section 2 recalls the general framework related to substructuring and the formulation of FETI method. After introducing the concept of error in constitutive relation, section 3 presents its formulation within domain decomposition. Then, section 4 details the specific implementation in the context of the FETI method. Finally, section 5 presents numerical assessments, especially to show that a good estimation can be obtained far earlier than the solver converged (in the sense of domain decomposition iterative solver).

FRAMEWORK AND BASICS

Reference mechanical problem

Let us consider the static equilibrium of a structure which occupies the open domain Ω ⊂ R d and which is submitted to given body forces f , to given traction forces g on ∂ f Ω and to given displacements u d on the complementary part ∂ u Ω ≠ ∅. We assume the structure undergoes small perturbations and that the material is linear elastic, characterized by the Hooke's tensor H.

For an open subset ω ⊂ Ω, and denoting by u the unknown displacement field, ε(u) the symmetric part of the gradient, σ the Cauchy stress tensor, we introduce the following subspaces Ka(Ω) of kinematically admissible and Sa(Ω) statically admissible fields:

Ka(ω) = v ∈ H 1 (ω) d such that tr(v) ∂uω = u d (1a) Sa(ω) = τ ∈ L 2 (ω) d 2 sym such that ω τ ∶ ε(v) dω = ω f.v dω + ∂ f ω g.v dΥ, ∀ v ∈ Ka 0 (ω) (1b)
where tr is the trace operator and KA 0 (ω) is the vector space associated to Ka(ω).

The mechanical problem set on Ω can then be formulated as:

Find (u, σ) ∈ Ca(Ω) × Sa(Ω) verifying: σ = H ∶ (u) (2)

Substructured formulation

Let us consider a decomposition of domain Ω in open subsets (Ω (s) ) 1⩽s⩽N sd (N sd is the number of subdomains) so that Ω (s) ∩ Ω (s ′ ) = ∅ for s ≠ s ′ and Ω = ∪ s Ω(s) . For any N sdturple (v (s) ) 1⩽s⩽N sd of local fields denoted by v ◻ , we define the global assembling operator A by:

v = A(v ◻ ) ⇔ v Ω (s) = v (s) (3)
Under this framework, kinematic and static admissibility on the whole structure can be restricted to each sub-structure Ω (s) provided the verification of interface conditions, namely displacements continuity and tractions balance (action-reaction principle):

A(u ◻ ) ∈ Ka(Ω) ⇔ u (s) ∈ Ka(Ω (s) ), ∀s tr(u (s) ) = tr(u (s ′ ) ) on Υ (ss ′ ) , ∀(s, s ′ ) (4a) A(σ ◻ ) ∈ Sa(Ω) ⇔ σ (s) ∈ Sa(Ω (s) ), ∀s σ (s) .n (s) + σ (s ′ ) .n (s ′ ) = 0 on Υ (ss ′ ) , ∀(s, s ′ ) (4b)
The sub-structured problem may then be formulated as:

Find (A(u ◻ ), A(σ ◻ )) ∈ Ca(Ω) × Sa(Ω) verifying: σ (s) = H ∶ u (s) (1 ≤ s ≤ N sd ) (5)

Finite element approximation

Let Ω h be a tessellation of Ω associated to a finite dimensional subspace Ka h (Ω) of Ka(Ω) and leading to a conforming substructuring so that (i) each element only belongs to one subdomain and (ii) nodes are matching on the interfaces. Under those assumptions, each degree of freedom is either located inside a subdomain (subscript i) or on its boundary Υ (s) = ∪ s ′ Υ (ss ′ ) (subscript b) where it is shared with at least one neighboring subdomain. Approximation u h ∈ Ka h (Ω) of u may then be expressed with the d × N dof matrix ϕ of shape functions (which form a basis of Ka h (Ω)) and the vector of nodal unknowns u as u h = ϕu.

Introducing λ (s) b , vector of unknown efforts imposed on the interface of a subdomain Ω (s) by its neighbors, the finite element problem can be written on each subdomain separating internal and boundary degrees of freedom (dof):

K (s) ii K (s) ib K (s) bi K (s) bb u (s) i u (s) b = f (s) i f (s) b + 0 (s) i λ (s) b (6)
where K (s) is the (symmetric positive definite) stiffness matrix of domain Ω (s)

h and f (s) is the vector of generalized forces. Eq. ( 6) then allow to express equilibrium on sub-domains in a condensed form function of the dof located on Υ (s) (eq. 7a). In addition to that, the discretization of interface conditions included in (4a) and (4b) can be expressed as dof-todof relationships thanks to the conforming feature of the decomposition. As a result, finite element approximation of the substructured reference problem leads to the following interface system:

S (s) u (s) b = b (s) + λ (s) b (7a) s A (s) u (s) b = 0 (7b) s A (s) λ (s) b = 0 (7c)
with:

S (s) = K (s) bb -K (s) bi K (s) ii -1 K (s) ib b (s) = f (s) b -K (s) bi K (s) ii -1 f (s) i
where A (s) and A (s) are assembling operators 1 so that A (s) enables to formulate the mechanical equilibrium of interfaces (4b) and A (s) enables to formulate the continuity of displacements (4a) (see [START_REF] Gosselet | Non-overlapping domain decomposition methods in structural mechanics[END_REF] for more an extensive description of all operators). Let's bring back the fundamental orthogonality property verified by assembling operators:

s A (s) A (s) T = 0 (8)

FETI approach

The starting point of the Finite Element Tearing and Interconnecting domain decomposition [START_REF] Farhat | Implicit parallel processing in structural mechanics[END_REF], is to introduce a unique interface traction unknown λ b so that interface equilibrium is always insured when expressing λ

(s) b : λ (s) b = A (s) T λ b ⇒ s A (s) λ (s) b = 0 (9)
Displacements can yet be deduced from λ b (10a) if it satisfies Fredholm's alternative (10b) on each substructure:

u (s) b = S (s) + b (s) + A (s) T λ b + R (s) b α (s) (10a) 0 = R (s) b T b (s) + A (s) T λ b (10b)
where α (s) is the unknown magnitude of rigid body motions and superscript + denotes the generalised inverse. The FETI solver then consists in iteratively finding an interface effort λ b (under the previous constraint) insuring the continuity of interface displacement

0 = ∑ s A (s) u (s)
b , which leads to the well-known interface problem:

S d G G T 0 λ b α = b d e (11)
where:

S d = N sd s=1 A (s) S +(s) A (s) T G = A (1) R (1) b . . . A (N sd ) R (N sd ) b α = α (1) T . . . α (N sd ) T T b d = - N sd s=1 A (s) S (s) + b (s) p e = -b (1) p
T R (1) b ⋅ ⋅ ⋅ -b (N sd ) p T R (N sd ) b T 1
In the case of two subdomains, one has

∑ s A (s) λ b = λ (1) b + λ (2) b = 0 and ∑ s A (s) u (s) b = u (1) b -u (2) b = 0.
In order to iteratively solve problem [START_REF] Le | Domain decomposition methods in computational mechanics[END_REF] in an efficient way, the best preconditioner S-1 d for operator S d is built as an assembling of inverses of local operators S (s) + previously involved:

S-1 d = N sd s=1 Ã(s) S (s) Ã(s) T (12)
where ( Ã(s) ) s are scaled assembling operators verifying the fundamental property (13),

whose typical example is the multiplicity scaling (14

) 2 : s A (s) Ã(s) T = I (13) Ã(s) T = A (s) T j A (j) A (j) T -1 (14)
Let us note that more elaborate operators taking the heterogeneity into account may be used in the case of heterogeneous structures [START_REF] Daniel | A simple and efficient extension of a class of substructure based preconditioners to heterogeneous structural mechanics problems[END_REF].

A POSTERIORI ERROR ESTIMATION

Error in constitutive relation

Let us come back to the continuous level and introduce the functional error in constutive relation measuring the non-verification of constitutive relation between a given pair of displacement and stress (v, τ ):

e rdc,ω (v, τ ) = τ -H ∶ ε(v) H -1 ,ω with • H -1 ,ω = ω • ∶ H -1 ∶ • dω 1 2 (15)
Recall that given the tessellation Ω h and its associated subspace of admissible displacements Ka h (Ω) ⊂ Ka(Ω), the usual finite element approximation of reference problem (2) writes:

Find u h ∈ Ca h (Ω) such that σ h = H ∶ ε(u h ) satisfy: Ω σ(u h ) ∶ ε(v h ) dΩ = Ω f.v h dΩ + ∂ f Ω g.v h dΥ, ∀ v h ∈ Ca 0 h (Ω) (16)
Then, on the basis of the solution (u h , σ h ) of ( 16) satisfying the classical properties u h ∈ Ka(Ω) and e cr,Ω (u h , σ h ) = 0 but σ h ∉ Sa(Ω), the basic principle of the error in constitutive relation consists in deducing an admissible displacement-stress pair (û h , σh ) ∈ Ka(Ω) × Sa(Ω) in order to measure the residual on the constitutive equation [START_REF] Parret-Fréaud | Fast estimation of discretization error for fe problems solved by domain decomposition[END_REF]. Hence, the evaluation of e cr,Ω (û h , σh ) ≥ 0 for any admissible pair (û h , σh ) provides a guaranteed upper bound of the global error thanks to the Prager-Synge theorem:

e H,Ω = u -u h H,Ω ≤ [e rdc,Ω (û h , σh )] 1 2 with • H,ω = ω ε(•) ∶ H ∶ ε(•) dω 1 2 (17)
2 In the case of two subdomains, one has ∑ s Ã(s) λ

(s) b = 1 2 λ (1) b -λ (2) b
.

Ka h (Ω) being a subspace of Ka(Ω), the construction of an admissible displacement field ûh is straightforward since it can be taken equal to u h . On the other hand, as σ h is not statically admissible, the construction of an admissible stress field σh ∈ Sa(Ω) is a crucial point which has already been widely studied in the literature. In addition to the use of dual formulation [START_REF] Kempeneers | Pure equilibrium tetrahedral finite elements for global error estimation by dual analysis[END_REF], classical methods consist in post-processing a statically admissible field from Field σ h (obtained by a displacement formulation) through the element equilibration techniques [START_REF] Ladeveze | Error estimate procedure in the finite element method and applications[END_REF][START_REF] Ladevèze | New advances on a posteriori error on constitutive relation in f.e. analysis[END_REF], which have been improved by the use of the concept of partition of unity [START_REF] Ladevèze | A new non-intrusive technique for the construction of admissible stress fields in model verification[END_REF] or flux-free-based methods [START_REF] Pares | Subdomain-based flux-free a posteriori error estimators[END_REF][START_REF] Gallimard | A constitutive relation error estimator based on traction-free recovery of the equilibrated stress[END_REF][START_REF] Moitinho De Almeida | Recovery of equilibrium on star patches using a partition of unity technique[END_REF]. In most cases they involve the computation of efforts on "star-patches" which are the set of elements sharing one node, for each node of the mesh. Though rather simple, these computations are in great number and thus expensive.

In the following, we note by F h the algorithm which has been chosen to build an admissible stress field σh . Whatever the choice, the algorithm takes as input not only the finite element stress field σ h but also the continuous representation of the imposed forces (f, g).

σh = F h (σ h , f, g) ∈ SaΩ
The algorithm we have used for our applications is the one proposed in [START_REF] Ladeveze | Error estimate procedure in the finite element method and applications[END_REF] using a three degrees higher polynomial basis when solving the local problems on elements [START_REF] Babuška | Validation of a posteriori error estimators by numerical approach[END_REF].

Formulation within a generic domain decomposition

Following the previously presented idea of error in constitutive relation in a sequential framework, the key point for an efficient evaluation in a substructured context (without overlapping) is to define admissible pairs (û

(s) h ,σ (s)
h ) ∈ Ka(Ω (s) ) × Sa(Ω (s) ) on each subdomain so that the associated assembled pair is admissible for the reference problem

(A(û ◻ h ), A(σ ◻ h )) ∈ Ka(Ω) × Sa(Ω).
Due to the absence of overlap, the additive structure of the associated error in constitutive relation leads to a fully parallel evaluation of the a posteriori error estimator:

e cr,Ω (A(û ◻ h ), A(σ ◻ h )) 2 = N sd s=1 e cr,Ω (s) (û (s) h , σ(s) h ) 2
However, the application of a classical recovery strategy to compute admissible fields raises two difficulties in a substructured context. First, the star-patches can not be employed on the boundary nodes without assuming communication between subdomains. Though these exchanges would remain limited, we propose an alternate strategy to achieve full parallelism without impairing the properties of the error in constitutive relation. Second, recall that the FETI method used on the discrete substructured problem (7) consists in iteratively solving eq. ( 11) so that displacement continuity (4a) is only verified (up to a certain precision) once the solver converged. Thus recovering strategies need to be adapted so that the local fields (û (s) h , σ(s) h ) satisfy the interface conditions.

PARALLEL RECOVERY OF ADMISSIBLE FIELDS

Kinematically admissible fields

In order to ensure interface Condition (4a) when building û(s) bh , the last condition can be directly transposed into a discrete counterpart thanks to the matching discretizations on the interface:

h ∈ Ka(Ω (s) ) so that A(û ◻ h ) ∈ Ka(Ω),
û(s) b = û(s ′ ) b , ∀(s, s ′ ) ⇐⇒ s A (s) û(s) b = 0 (18)
However, displacement continuity (7b) is only verified upon convergence in the FETI method, but local displacements on each substructure may be computed at each iteration through (10a) since condition (10b) is always satisfied. The associated gap of interface displacement is given by:

u b = s A (s) u (s) b (19)
Then, we choose to use the scaled assembly operators introduced for the preconditionning step [START_REF] Mandel | Balancing domain decomposition[END_REF], and define û(s) b as:

û(s) b = u (s) b - Ã(s) T u b (20)
Thereafter, the property (13) satisfied by operators Ã(s) ensures the verification of the discrete continuity condition (18) on û(s) b . Finally, one finite element problem solved independently on each subdomain with imposed Dirichlet conditions on the interface given by û(s) b enables to deduce the kinematically admissible displacement field û = A(û ◻ h ) ∈ Ka(Ω):

û(s) i = K (s) ii -1 f (s) i -K (s) ib A (s) T û(s) b (21) û(s) h = ϕ (s) û(s) = ϕ (s) i ϕ (s) b û(s) i û(s) b
Let us note a key feature of this approach which comes from the fact that all the involved operations are already realized during the steps related to the Dirichlet's preconditioner [START_REF] Mandel | Balancing domain decomposition[END_REF], so that all finite element quantities (even the internal ones) are available at no cost: the quantity u b is directly available during the classical solution procedure (without computing any α (j) ) which is based on an initialization/projection algorithm [START_REF] Farhat | Implicit parallel processing in structural mechanics[END_REF], and the displacement field u (s) can be defined up to an element of the kernel (a rigid body motion) since only its symmetric gradient is used during the computation of the error.

A. Parret-Freaud, P. Gosselet and C. Rey

Statically admissible fields

In order to ensure interface Condition (4b) when building σ(s) h ∈ Sa(Ω (s) ) so that A(σ ◻ h ) ∈ Sa(Ω), a generic way consists in introducing on each subdomain a continuous balanced interface traction field F (s) bh defined on Υ (s) which satisfy:

F (s) bh + F (s ′ ) bh = 0 sur Γ (ss ′ ) (22a)
Ω (s) f (s) .ρdΩ + ∂ f Ω (s)
g (s) .ρdS +

Γ (s) F (s) bh .ρdS = 0 ∀ρ ∈ KaR(Ω (s) ) (22b)
where KaR(Ω (s) ) is the set of rigid body motions which are compatible with Dirichlet conditions imposed on ∂ u Ω (s) , and (22b) stands for the translation of Fredholm's alternative.

We then associate F (s) bh to a finite element nodal reaction field

λ(s) b through an algorithm F (s) bh = G h λ(s) b
requiring that the discrete field λ(s) b and the continuous field F (s) bh develop the same virtual work in any finite element displacement field:

Υ (ss ′ ) F (s) bh ⋅ ϕ (s) j Υ (ss ′ ) dS = λ(s) b,j (23)
where j denotes a node of the interface, ϕ (s) j its associated shape function and λ(s) b,j the corresponding nodal component of λ(s) b . We get the following discrete counterparts for λ(s) b :

s A (s) λ(s) b = 0 R (s) b T λ(s) b + b (s) = 0 (24)
In the FETI solver, the nodal interface fields are by construction ( 9) always balanced at the interface and associated to well-posed discrete Neumann problems on each substructure thanks to (10b), allowing us to directly set λ(s

) b = λ (s) b . Setting up the algorithm G h in order to get F (s) bh from λ(s) b
is strongly linked to a choice of interpolation basis. A simple choice is then to use the basis of finite element shape functions on Υ (s) :

Fbh = ϕ (s) Υ (s) F(s)
b , which, though leading to a coarse interpolation, doesn't alter the global performance of the method 3 .

Once the traction field F (s) bh has been recovered, we get the statically admissible field σh = A(σ ◻ h ) ∈ Sa(Ω) through a two-step process. First one is to solve a finite element problem (25a) on each subdomain with imposed Neumann conditions on the interface given by F (s) bh . The next one is to apply a sequential recovery algorithm F h on each finite element stress field H ∶ ε(ϕ (s) ũ(s) ) previously obtained (25a) taking into account F (s) bh as an imposed force on Υ (s) :

ũ(s) = K (s) + f (s) + t (s) T λ(s) b (25a) σ(s) h = F eq H ∶ ε(ϕ (s) ũ(s) ), f (s) , g (s) , G h ( λ(s) b )) (25b) (25c)
It has to be observed that the fully parallel procedure G h proposed above leads to a different admissible traction field than the one which would have been obtained using standard patch-technique [START_REF] Ladevèze | New advances on a posteriori error on constitutive relation in f.e. analysis[END_REF] (referred in the sequel as the sequential approach) so that the parallel error estimation is different from the standard sequential one. Finally, let us note that the only extra operations are due to the use of algorithms G h (to compute Fbh ) and F h (to compute σh ).

NUMERICAL ASSESSMENT

The performances of our parallel error estimator has been assessed on a 2D toy problem of a Γ-shape structure clamped on its basis and submitted to traction and shear on its upper-right side, with an isotropic linear elastic material behavior. A set of computations with increasing number of subdomains has been driven on a sequence of regular meshes generated with a increasing refinement. In addition to the new parallel error estimator e ddm cr , the standard sequential one e seq cr and the true error e h (using a reference field u ex on a very fine mesh) has been computed:

e seq cr = e cr,Ω (û h , σh ) e ddm cr = e cr,Ω (A(û ◻ h ), A(σ ◻ h )) e h = u ex -ûh H,Ω = u ex 2 H,Ω -ûh 2 H,Ω (26)
At last, let us note that all our computations are driven in the ZeBuLoN finite element code [START_REF]Z-set user manual[END_REF], using elements of polynomial degree p = 1. The results when convergence of the domain decomposition solver is reached, which can be found in the paper [START_REF] Parret-Fréaud | Fast estimation of discretization error for fe problems solved by domain decomposition[END_REF], show that e ddm cr barely depends on the substructuring so that the estimates are quite similar whether they are conducted on a single domain or on N sd subdomains. As a conclusion, the parallel error estimator e ddm cr enables to recover the same efficiency factor as the standard sequential one, while the CPU-time is divided by N sd .

Another interesting feature associated to the use of an iterative solver for the domain decomposition (DD) problem is that the discretization error estimation can be conducted before DD convergence is reached, that is in presence of displacement discontinuity at the interface when using FETI. We illustrate this point on a small crack opening problem (fig. 1) already used in other papers [START_REF] Pares | Subdomain-based flux-free a posteriori error estimators[END_REF][START_REF] Ladevèze | A new non-intrusive technique for the construction of admissible stress fields in model verification[END_REF]. The behavior is supposed to be isotropic linear elastic (E = 2000 Pa and µ = 0.3) and plane stress hypothesis is assumed. Beside, the FETI algorithm used on the 14-sub-domain decomposition is equipped with the Dirichlet preconditioner (12) while its convergence criterion (which stands here for the interface displacement gap) is set to 10 -6 . Figure 2 shows the convergence curve of e ddm cr during the FETI iterations, plotted as a function of the corresponding residual that is the normalized displacement gap at the interface. The curve show a rapid convergence of the parallel error estimator along iterations of the solver, so that e ddm cr can be considered as converged when the residual reaches an order of magnitude of 5.10 -2 which corresponds to at most 4 iterations whereas the solver convergence is achieved in about 15 iterations. The "L"-shaped curves show that the impact of residual of the DD solver is preponderant only at the first iterations (when interface fields are very poorly estimated), afterwards e ddm cr stagnates at a value close to e seq cr which is only associated to the discretization error. The maps of the elementary contributions e ddm cr,E to the parallel error estimator e ddm cr at iterations 1, 2 and 3 corroborate the previous observation. An interesting fallout of the foregoing is that when willing to carry out remeshing procedures, the maps obtained after few iterations of the solver are sufficient to define correct refinement instructions.

CONCLUSIONS

In this paper, we presented a new method to handle error estimation within domain decomposition approaches. Our approach relies on the construction of admissible interface fields followed by a fully parallel construction of a displacement-stress pair that is kinematically and statically admissible on the whole structure. We then obtain a simple and efficient extension to error in constitutive relation in a substructured context. Moreover, this construction can be realized at every iteration of the domain decomposition solver, at a reduced numerical cost when taking advantage of the associated preconditionners. As a result, not only our method enables to divide the huge CPU cost associated to error estimation, but a satisfactory estimate can also be accessed far before the solver reaches convergence. Indeed, results show that even roughly estimated interface fields enable to obtain a good estimation of the discretization error and correct maps of elementary contributions. This last feature let us envisage, on the one hand, the development of a new convergence criterion for the iterative solver, and, on the other hand, the set up of efficient mesh adaptation procedures relying on estimates obtained after only few iterations. Such procedures are currently studied and leads to different approaches whether the remeshing step may preserve or not the topology of the sub-structuring. Other current studies deals with heterogeneous problems, whose results will be presented during the conference.

  (a) Finite-element problem (b) Substructuring (N sd = 14) (c) Elementary contributions e ddm cr,E at convergence

Figure 1 :

 1 Figure 1: Crack opening problem

Figure 2 :

 2 Figure 2: Convergence of e ddm cr vs. interface residual and elementary contributions to e ddm cr

Let us note that more consistent choices needs to allow discontinuities at nodes, a counterpart of which is to use an additional minimization problem. Some alternative choices are currently studied and will be discussed.