J Comte

S Legand

Determination of carbon-14 source term in French SIERs

Carbon Organic mass balance is calculated • ≈ 12 % is identified in LiOH "blank", virgin resins and SIERs solutions.

• Only one organic compound identified in resin desorption solutions at very low concentration ≈ 1 % of Total Organic Carbon (TOC).

• Identified carboxylic acids do not complete the carbon mass balance, even in LiOH solution. So, other molecules are in solution.

, Contamination in LiOH makes difficult the determination of organic molecules in resin desorption solutions.

Introduction / Context

Task 4-2 : Aim of the project achieved Two analytical methodologies were developed for total 14 C and speciation measurement and applied on real samples (Dry/Wet) Ion Exchange Resins (IERs as combination of cationic, anionic and/or mixed bed resins) are used for purification of the various circuits in various operating history.

Results show

• EDF-SIERs are classified as Low and Intermediate level Short Lived waste (LILW-SL)

• Contain short-lived radionuclides, such as cobalt-60, and may also contain strictly limited amounts of long-lived radionuclides such as carbon 14.

• Are disposed off at the Andra CSFMA surface waste disposal facility after packaging

Spent ion exchange resins samples

• IER (0.2 to 0.

EDF-B sample : inconsistent results

Cutting of a SIERs package

  That it is difficult to compare different wet SIERs (H 2 O content variability / analytical problems) -14 C speciation on one sample showing ~75 % of mineral carbon But question remain about the representativeness of the results => Origin of spent ion exchange resins EDF operates 58 PWR in France. Spent Ion Exchange Resins (SIERs) are discharged under water and stored before final conditioning. Discharge conditions criteria are usually chemical or radiochemical saturation, maximum service lifetime... Determination of 14 C source term (WP4 Task 4

  5 g) mixed with carrier (graphite) • Combustion under O 2 (25 bars) • CO 2 produced by mineralization : -trapped in NaOH solution as CO 3 2-in the mineralization device during combustion and cooling -to achieve the recovery of 14 C, the mineralization system is degassed in a venting station and the gases are flushed into NaOH bubblers Sample + carbon carrier Sodium hydroxyde solution (Inside) • IER (0.2 to 0.5 g) spiked by carbonate solution • First step : measurement of mineral 14 C -Acidification by H 2 SO 4 -Degassing under N 2 and vacuum into NaOH bubblers -Measurements of CO 3 2-carrier (recovery) and 14 C • Second step : measurement of organic 14 C -Concentrated acid (H 2 SO 4 /HNO 3 ) with oxidant (H 2 O 2 ) -Degassing under N 2 and vacuum reflux into NaOH bubblers -Measurement 14 C Similar as « Magnusson method» thesis 2007 Optimizations performed on virgin resins spiked with 14 C as carbonate and/or glucose Combustion with

  Samples are representatives of the waste before conditioning (not of the "packaged waste") Samples are made of SIERs from different circuits depending on the sampling campaignSampling of SIERs (~10g) are performed during this process and transferred to CEA laboratorySIERs packaging "Mercure" process conditioningMobile unit on EDF site Desorption experiments (WP4 Task 4-3) Desorption experiments performed to identify organic compounds fixed on the SIERs • Experimental methodology -Desorption~1 g of wet SIERs by 10 ml of LiOH 10 -3 up to 1 Mol/l • Carbon 14 results -~ 10 % of 14 C remains in SIERs after desorption experiment (in organic form) -Desorption 14 C balance (~ 1400 Bq/g) does not fit with total carbon 14 (2125 Bq/g with the combustion method)  carbon 14 loss in the gas phase ? This work was carried out under the CAST project : CArbon 14 Source Term Work package 4 with partial funding from the European Union's Seventh Framework and from the programme for research, technological, development and demonstration under grant agreement no. 604779, ANDRA, EDF and the CEA. The authors would also like to thank EDF for having provided the ion exchange resin samples.  problem is identified in the case of drying => more experiments are needed  Samples representativeness is questioned => a statistical approach is needed to take into account the origin of the SIERs (various campaigns) Task 4-3 : Aim of the project achieved Desorption experiments permit to obtain aqueous solutions containing 14 C in organic form for the determination of the hydro soluble organic compounds Measurements with different techniques were used Results show One organic compound representative of resin desorption is observed but doesn't complete the total carbon organic mass balance ≈ 1 % of TOC identified This poster presents the methodologies and results of experiments carried out on the spent ion exchange resins (SIERS) from EDF PWR (pressurized water reactor) at the Commissariat à l'Énergie Atomique et aux Énergies Alternatives (CEA beneficiary N°4 ) within the framework of the European program Carbon-14 Source Term (CAST) Work Package 4. The studies concern both carbon-14 measurements and its speciation (inventory) and the measurement of organic molecules issued from extraction experiments on SIERS. Two analytical methodologies are needed to check the accuracy of measurements Application on drying SIERs show that some inconsistent results were obtained Carboxylic acids by anionic chromatography Low molecular weight volatile molecules by GC-MS Amines by cationic chromatography Molecules with higher molecular mass by ESI-MS Acetate & formiate detected in LiOH solutions (contamination), only glycolate come from resin desorption (< 8 mg acid .L -1 ) Derivatization reaction performed to improve volatile organic compounds identification, nothing detected on resin desorption No amine detected with the developed method, no chemical degradation of functional groups observed in these conditions LiOH matrix but no organic compound with higher molecular weight in solution is detected (probably salts hide signal) Interpretation of organic carbon speciation in resin desorption solutions Analytical development and sample analysis for the determination of the hydro soluble organic compounds * GC-MS : Gas Chromatography -Mass Spectrometry -ESI-MS -: Electrospray -Mass Spectrometry Total carbon determined by combustion -Speciation determined by acidic dissolution

Name of the sample Power reactor unit (MWe) Type of fuel Sampling date Received date

						Mass (g)	Dose rate
	EDF-A	1450	UOX	01/2014	06/2014	10 15 µSv/h (10g)
	EDF-B	1300	UOX	11/2014	02/2015	4

*10 300 µSv/h (10g)


Chao Zheng

Oleg Rubel
email: rubelo@mcmaster.ca

Mikael Kepenekian

Xavier Rocquefelte

Claudine Katan
email: claudine.katan@univ-rennes1.fr

Mikaël Kepenekian

Electronic properties of Pb-I deficient lead halide perovskites

teaching and research institutions in France or abroad, or from public or private research centers.

I. INTRODUCTION

Recently, a new family of hybrid-halide perovskites called "deficient" perovskites (dperovskites) and "hollow" perovskites was reported [START_REF] Leblanc | Leadand iodide-deficient (CH 3 NH 3 )PbI 3 (d-MAPI): The bridge between 2D and 3D hybrid perovskites[END_REF][START_REF] Ke | Enhanced photovoltaic performance and stability with a new type of hollow 3d perovskite {en}FASnI 3[END_REF][START_REF] Spanopoulos | Unraveling the chemical nature of the 3D "hollow" hybrid halide perovskites[END_REF][START_REF] Ke | Efficient lead-free solar cells based on hollow {en}MASnI 3 perovskites[END_REF][START_REF] Ke | Diammonium cations in the FASnI 3 perovskite structure lead to lower dark currents and more efficient solar cells[END_REF][6] . Their structure is different from conventional AM X 3 perovskites where a relatively small organic cation A + , such as methyl ammonium (A + = MA + ), is centered within an inorganic M X - 3 cage formed by lead or tin halides. Deficient perovskites have missing (M -X) + units with corresponding voids filled by larger cations, such as hydroxyethylammonium (A + = HEA + ), ethylenediammonium (A + = EDA + ), or thioethylammonium (A + = TEA + ). The defect (d) can be expressed using Kröger-Vink notations as

d × = v Pb + v • I + (A ) • Pb-I . (1)
The details of the defect reaction in Eq. ( 1) indicates that d consists of two intrinsic defects and one extrinsic defect. The superscript × indicates null charge for the current site. v Pb is the first intrinsic defect which signifies the vacancy of Pb. The double prime denotes a negative charge (-2) due to missing of a Pb 2+ . The second intrinsic defect is the vacancy of I v • I , where the superscript • denotes a net single positive charge (+1) due to missing of I -. (A ) • Pb-I is the extrinsic defect which demonstrates that A + cation occupies the vacancies of Pb and I. This defect contributes a net single positive charge due to the presence of A + .

The general formula of d-perovskites is (A, A ) 1+x M 1-x X 3-x 1 . The concentration x of such structural defects can be quite large (the range of 0 < x < 0.2 has been demonstrated

experimentally [START_REF] Leblanc | Leadand iodide-deficient (CH 3 NH 3 )PbI 3 (d-MAPI): The bridge between 2D and 3D hybrid perovskites[END_REF] ) leading to a range of 3 < X/M < 3.5 that captures a deviation from conventional stoichiometry of perovskites. The difference between deficient and hollow perovskites comes from an ordered vs stochastic arrangement of defects, respectively. Despite a large concentration of structural defects, the material retains its 3D architecture and can be considered as a transition between the 2D and 3D structures [START_REF] Katan | Quantum and dielectric confinement effects in lowerdimensional hybrid perovskite semiconductors[END_REF] .

The most noticeable change in electronic properties of d-perovskites is the opening of an optical band gap. The band gap of d-MAPbI 3 with I/Pb ratio of 3.5 (x = 0.2) increases from 1.55 to 2.15 eV 1 . Spanopoulos et al. [START_REF] Spanopoulos | Unraveling the chemical nature of the 3D "hollow" hybrid halide perovskites[END_REF] also reported increase of the band gap from 1.52 to 1.85 eV for hollow MAPbI 3 perovskites (I/Pb ratio = 3.6) and from 1.25 to 1.42 eV for hollow MASnI 3 perovskites (I/Sn ratio = 3.6). Discrepancies in the experimental data can be attributed to different large cations (HEA + vs EDA + ). One should also note that the d-perovskites have been reported to present a dual band gap marked by the existence of two shoulders in the UV/vis spectra [see Fig. 2(a) in Ref. 1] that leads to another uncertainty. Other noticeable properties of deficient/hollow perovskites are improved stability [START_REF] Ke | Enhanced photovoltaic performance and stability with a new type of hollow 3d perovskite {en}FASnI 3[END_REF][START_REF] Spanopoulos | Unraveling the chemical nature of the 3D "hollow" hybrid halide perovskites[END_REF]6 , improved photoluminescence efficiency of hollow MASnI

II. METHOD

Density functional theory 12,13 (DFT) calculations were performed using the Vienna ab initio simulation package 14 (VASP) and projector augmented-wave potentials [14][15][16] The Pb I framework was kept fixed at experimental positions and lattice parameters 1 .

This implies that the same lattice constant and atomic distortions apply to all compositions, including the stoichiometric phase of MAPbI 3 [see supporting information (SI), Table SI].

The experiment 1 indicates only a small increase of the lattice parameters for d-MAPbI 3 .

Cesium atoms were used instead of MA and HEA cations (Fig. 2) to reduce computational efforts. This approach is justified as long as we preserve experimental positions for the leadiodide framework, which determines the band gap and dispersion near to the band edges. In order to ensure that the conclusions drawn in this paper are not artifacts of the approximation made regarding the structure of d-perovskites, we investigate a structure where HEA + organic cation is substituted in place of (Pb I) + unit in tetragonal MAPbI in the range 0-1 for each n's energy eigenstate to reflect its Bloch k character. The number of spectral weights computed per eigenstate corresponds to the multiplicity of the supercell.

The approach to the calculation of spectral weights is based on remapping the supercell reciprocal space with a mesh that is compatible with the translational symmetry of a primitive cell. Computed Bloch spectral weights for each eigenstate in a supercell calculation fulfill the normalization k w n (k) = 1, which implies that the cumulative probability of finding an electron in any of multiple k's states adds up to one. Further details of the method and implementation can be found elsewhere 20,21 .

Vesta 22 was used for visualization of atomic structures. Structure files used in this work can be accessed at the Cambridge Crystallographic Data Center (CCDC) under deposition

No . 1962341-1962347 and 1962349.

III. RESULTS AND DISCUSSION

A. Structure of d-perovskites and defect formation enthalpy

It is worth commenting on the structure of d-perovskites before going into discussion of their electronic properties. Here we adopt the structure proposed in Ref. 1, where the lattice parameters and coordinates of Pb I framework are specified. Fig. 2 shows the stoichiometric structure as well as two d-perovskite structures that differ by the I/Pb ratio (3.22 and 3.5) corresponding to the defect concentration of x = 0.1 and 0.2, respectively. The closer the ratio to its stoichiometric value of 3, the less deficient is the perovskite structure. The ratio of I/Pb = 3.5 is at the upper limit of deficiency attained experimentally 1 , if not hollow. In the sections that follows, Cs atoms are used instead of organic cations (MA or HEA) for computational performance reasons (see also the case study with HEA in Method section).

The self-compensation defect formed by A cation substituted for M -X unit might also be unintentionally present in smaller concentrations in stoichiometric AM X 3 perovskite structures. To evaluate a formation energy of the corresponding defect [Eq. ( 1)] we built a (2)

Giving the size of the supercell, it corresponds to a defect concentration of x = 0.125.

We assign the DFT total energy difference between reactants and products to the defect formation enthalpy

∆H d = E MA 8 Pb 8 I 24 + E PbI 2 -E MA 8 (HEA)Pb 7 I 23 -E HEAI . (3)
The calculations yield ∆H d = 0.78 eV. This result can be compared with another common Schottky defect in perovskites, namely PbI 2 vacancy. The later has a formation energy of about 0.2 eV 23 suggesting that the defect described by Eq. ( 1) is less likely to form spontaneously.

B. Band gap opening, band alignment, dual band gap and optical properties

There is a significant opening of the band gap in the deficient structures relative to the reference (stoichiometric) perovskite. Computed band gaps are listed in Table I. The band gap opening is more prominent when relativistic effects are included. The calculated band gap increment (0.74 eV with SOC and 0.36 eV without SOC) is within the range of the experimental band gap opening quoted in Sec. I.

It is quite unusual that defects lead to opening of the band gap. It implies that the defect described in Eq. ( 1) does not create states within the band gap. The Pb I vacancy introduces five undercoordinated I atoms and one Pb atom. A simplified orbital energy diagram of halide perovskites in Fig. 3(a) is intended to rationalize this result. Since both top of the valence band and bottom of the conduction band are formed by antibonding orbitals 24,25 , compensated dangling bonds of I -are expected to appear within the bulk of valence band states as confirmed by the analysis of atom-specific contribution to energy eigenstates (see SI, Fig. S3).

Opening of the band gap can be viewed as a confinement effect. Defects form local regions of A X with a wide band gap and confine electronic states within the perovskite region as shown in Fig. 3(b,c). The higher the concentration of defects, the smaller the distance between A X regions leading to a higher confinement energy. Since it is a confinement 6 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5127513 effect, we would expect it to be operational at a high enough concentration of defects when the distance between them becomes comparable to the lattice spacing.

It is interesting to see whether this type of defects will act as trap for charge carriers.

The greater band gap of d-perovskites alone does not provide an answer to this question, as there are different types of band alignment possible (namely type I or II, Fig. 1). The can be interpreted as the second gap. The monotonic behavior could be attributed to the growing density of states as the energy moves away from band edges into the bands (see SI, Fig. S4). Thus, to rationalize experimental observation, one needs to assume a coexistence of stoichiometric MAPbI 3 regions that contribute to E

(1) g

and Pb I deficient regions with higher optical transition energies at E

g as shown in Fig. 4(a).

C. Carrier transport properties

Since d-perovskites are considered as a potential absorber material in solar cells, [START_REF] Leblanc | Leadand iodide-deficient (CH 3 NH 3 )PbI 3 (d-MAPI): The bridge between 2D and 3D hybrid perovskites[END_REF][START_REF] Ke | Diammonium cations in the FASnI 3 perovskite structure lead to lower dark currents and more efficient solar cells[END_REF] it is interesting to explore how deficiencies affect intrinsic transport characteristics of this material. The band dispersion near CBE and VBE is shown in Fig. 6. As a general trend, d-perovskites show less dispersive bands, which was also noticed previously This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS

D. Band unfolding of d-perovskite

To gain further insight into changes of the band structure introduced by Pb I deficiency we employ a band unfolding technique. It requires the supercell to have lattice vectors aligned with lattice vectors of the primitive cell. Here we use the quasi-cubic cell (blue area on Fig. 7) as a primitive unit for projecting an unfolded band structure. We created a supercell in such a way that its lattice vectors are aligned with those of the primitive cubic cell while preserving the periodicity of the deficient perovskite structure at the same time (Fig. 7). The smallest supercell that fulfills these conditions is shown as a pink area in Bloch character near R point implies that states are not localized despite the presence of structural defects. This is further corroborated when using inverse participation ratio 38 of electronic eigenstates as a measure of localization (see SI for details). In fact, analysis of the inverse participation ratio (see SI, Fig. S5) did not reveal any localized states near to the band edges either. This unusual behavior of d-perovskites is another manifestation of a defect tolerance inherent to hybrid halide perovskites 10,39 .

IV. CONCLUSION

The electronic structure evolution of d-perovskites with a general formula (A, A ) 1+x M 1-x X 3-x was investigated using the density functional theory for an intermediate and a high concentration of defects (x = 0.1 and 0.2, respectively). The formation enthalpy of Pb I deficient defects filled with a larger HEA + cation for charge compensation is about 0.8 eV, which is relatively large when compared to other native defects in halide perovskites. Our calculations reproduce the opening of a band gap also observed in previous theoretical This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. behave as a 3D material with anisotropic effective masses. As the density of defects increases (I/Pb = 3.5, x = 0.2), hollow regions connect into tunnels along the z axis, which remains the only direction for efficient band dispersion leading to a 1D electronic structure. Also, the band gap opening of d-perovskites exceeds that of 2D perovskites. This statement builds on the fact that a PL emission wavelength of 2D 1 ML MAPbI 3 perovskite is centered at 720 nm 40 , while d-MAPbI 3 have the secondary gap at shorter wavelengths (about 600 nm 1 ).

PLEASE CITE THIS ARTICLE AS

This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5127513

This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5127513

This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. FIG. 1: Band offsets at a semiconductor heterojunction illustrated for two types of band alignment . In type I, the conduction band of the well is lower than that of the barrier.

PLEASE CITE THIS ARTICLE AS

Oppositely, the valence band of the well lies higher than that of the barrier. In type II, both the conduction and valence bands of one of the two materials building the heterojunction are located below those of the second material.

3 4 ,

 4 and improved performance of Snbased perovskite solar cells[START_REF] Ke | Diammonium cations in the FASnI 3 perovskite structure lead to lower dark currents and more efficient solar cells[END_REF] .If the issue of defects in halide perovskites and their impact on photovoltaic or light emission performances has been intensively investigated (see for instance Refs. 8,9 and references herein), theoretical studies of d-perovskites are scarce and limited to Refs.3,5,6. Ab initio electronic structure calculations corroborate experimentally observed opening of the band gap; the band structure remains direct even at a high concentration of defects. It was also noticed that bands become less dispersed in hollow and d-perovskites[START_REF] Spanopoulos | Unraveling the chemical nature of the 3D "hollow" hybrid halide perovskites[END_REF]6 . The valence band dispersion in a lateral direction (i.e., perpendicular to defect tunnels) practically vanished in d-FAPbI 3 [Fig. 5(b) in Ref. 6] at a high concentration of defects (x = 0.2) that corresponds to formation of a continuous chain of defects (tunnels). Leblanc et al. 6 also pointed that the band gap opening with increasing x occurs mostly due to a shift of the conduction band edge (CBE) upward in energy, while the valence band edge (VBE) shifts down by approximately 0.1 eV. Similar to semiconductor heterojunctions induced band offset, defects will introduce energy misalignment. Fig. 1 illustrated two different types of band offsets, namely the type I band alignment or straddling gap and the type II band alignment or staggered gap. Identification of the type of band alignment is critical for investigating carrier transport properties in semiconductors. Here, we use first-principle calculations to investigate the electronic structure evolution of d-MAPbI 3 . Two deficient structures are selected with the ratio I/Pb = 3.22 and 3.5 to represent moderate (x = 0.1) and high density of defects (x = 0.2), respectively. We quantify changes in the band gap, band alignment, effective mass, and effective band structure of d-MAPbI 3 based on the structural model from Ref. 1. We anticipate the results to provide further insight into understanding a defect tolerance of perovskites 10,11 .

2 × 2 ×

 2 2 supercell of d-MAPbI 3 with and without defect as shown in Fig. S1 (see SI). The defect structure corresponds to the larger HEA + cation substituted in place of a Pb I unit This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5127513 of MAPbI 3 . The corresponding reaction can be expressed as MA 8 Pb 8 I 24 + HEAI → MA 8 (HEA)Pb 7 I 23 + PbI 2 .

  band alignment calculations were performed using core Cs-5s states as an energy reference to relate position of the valence and conduction bands in different calculations.(See Refs. 26,27 for details of this method as well as discussion of its accuracy and validation.) Our calculations suggest an energy misalignment of approximately δ v = 0.1 eV (TableI) between the VBE of stoichiometric and d-perovskite with VBE of d-perovskite being lower in energy as shown in Fig.4(a). This result is consistent with the recent data on d-FAPbI 3 6 . Most of the band gap difference between stoichiometric and d-perovskite comes from the CBE offset δ c (TableI). This result indicates a type-I band alignment between stoichiometric and d-perovskites. Thus, Pb I deficient regions should pose no threat for trapping of charge carriers. It is an interesting result, since defects are generally expected to create either shallow or deep traps within the band gap, which is not the case here.The disparity between the CBE and VBE band offset energies (δ c and δ v ) can be attributed to a much greater confinement energy for electrons than holes in the A X/AM X 3 heterostructure shown in Fig.3(b,c). The diagram was constructed based on experimental data for an electron affinity of 0.63 eV 28 and an ionization energy of 7.1 eV 29 for CsI as well as an electron affinity of 4.4 eV and an ionization energy of 5.9 eV for MAPbI 3 25 . The band diagram in Fig. 4(a) can now be used to explain an apparent dual band gap in reflectance measurements [see Fig. 2(a) in Ref. 1 also shown schematically in Fig. 4(b)]. The dual band gap character becomes more prominent with increasing deficiency of the perovskite. The first and second experimental gaps amount to E the possible origin of the dual band gap, we compute the optical response of d-perovskites. Calculations of optical properties (Fig. 5) indicate an anisotropy of the dielectric function in d-perovskites. The anisotropy becomes more prominent as the deficiency increases. The onset of absorption in d-perovskites is shifted to higher energies following the band gap trend and in accord with experimental observations 1,3 . The calculated imaginary part of the dielectric function (Fig. 5) shows a clear blue shift, but no features that This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5127513

  [START_REF] Spanopoulos | Unraveling the chemical nature of the 3D "hollow" hybrid halide perovskites[END_REF] . The band dispersion is mostly affected within the x-y plane. It is the plane perpendicular to direction of the tunnel formed after removing the Pb I unit [Fig.2(c)]. At high concentration of defects, the dispersion remains only along z axis that electronically resembles 1D structures.Effective mass calculations (TableII) corroborate conclusions drawn from the band dispersion data in Fig.6. The effective mass magnitude |m * /m e | in the stoichiometric perovskite is of the order of 0.09 for electrons and holes. Our calculated effective masses are about twice less that previous calculations for CsPbI 3 , |m * /m e | = 0.15 -0.20 30,31 and for MAPbI 3 |m * /m e | = 0.13-0.22 32 , and are also underestimated when compared to more accurate results based on many body perturbation theory within the GW approximation 33,34 . It is expected that the bare DFT-PBE will underestimates effective masses 35 due to the band gap error. Meanwhile, our computed reduced effective masses (0.08 -0.09) compare quite well to the experimental exciton reduced mass of MAPbI 3 0.104|m e | 36 . In d-perovskites, the effective masses of electrons and holes increases as the material becomes more defective. At the intermediate concentration of defects (the I/Pb ratio of 3.22, x = 0.1), the increase is moderate (by a factor of 4 at most). The effective mass becomes anisotropic. The lightest effective mass corresponds to the band dispersion along channels. The charge transport (especially holes) becomes effectively 1D when the concentration of defects (Pb I vacancies) approached the upper limit (I/Pb = 3.5, x = 0.2) and defects form a hollow tunnel along z axis.

  DOI:10.1063/1.5127513

Fig. 7 .

 7 Fig. 7. Its multiplicity factor is 5 × 5 × 2 relative to the cubic primitive cell. The unfolded band structures are shown in Fig. 8. A band structure of the stoichiometric phase was also unfolded for comparison. It shows a direct band gap at R point [Fig. 8(a)] as expected from previous studies of 3D lead-halide perovskites 24,37 . The Bloch character of electronic states near to the band edges of the stoichiometric perovskite is almost 100%. In the d-perovskite with high density of defects (I/Pb = 3.5, x = 0.2), 70 -80% of the Bloch character is preserved [Fig. 8(b)]. The band structure remains direct at R point. The band dispersion characteristics reflect the above discussion of effective masses. A well-preserved

  DOI:10.1063/1.5127513studies. The band gap opening is achieved mostly due to moving the conduction band edge upward in energy. In addition, the valence band edge shifts slightly toward lower energies by approximately 0.1 eV. This leads to the formation of a type-I band alignment between d-perovskites and the stoichiometric phase. As a result, Pb I deficient regions do not act as traps for charge carriers. This result is rationalized via a confinement model where the defect is viewed locally as a wide band gap region of A X.Presence of this type of defects in moderate quantities (x ≤ 0.1) should not hinder the ability of perovskite material to function as a good solar cell absorber. Based on effective mass calculations, we expect a reduced mobility of electrons and holes in d-MAPbI 3 (x = 0.1) by a factor of 4 at most. However, at high concentration of defects (x = 0.2), the charge transport ceases within a lateral plane due to an extremely heavy effective mass of holes thereby posing an issue for extracting photogenerated charge carriers. Interestingly, no localized states were detected near to the band edges of d-perovskites. It is attributed to the fact that both conduction and valence band edges are derived from antibonding orbitals.Investigation of the optical properties of d-perovskites reveal no anomalies in absorption spectra. The absorption spectrum increases monotonously up to a photon energy of about 2.5 eV. Thus, dual band gap features observed experimentally at a high concentration of defects (x > 0.15) are attributed to coexistence of deficient and stoichiometric regions.

Finally, we would

  like to comment on the electronic structure dimensionality of dperovskites. Deficient perovskites with a moderate density of defects (I/Pb = 3.22, x = 0.1)

FIG. 2 :FIG. 3 :FIG. 4 :

 234 FIG. 2: (a) Stoichiometric perovskite structure without deficiencies, (b) deficient perovskite with I/Pb ratio of 3.22 (x = 0.1), (c) I/Pb ratio of 3.5 (x = 0.2). Cesium atoms are used instead of MA cations. Defects are arranged in continuous channels along c axis.

FIG. 5 : 2 .FIG. 6 :FIG. 7 :FIG. 8 :

 52678 FIG. 5: Frequency-dependent imaginary part of the dielectric function for a stoichiometric perovskite structure without deficiencies (x = 0, violet) and d-perovskites with I/Pb ratio of 3.22 and 3.5 (x = 0.1, green, and 0.2, orange, respectively). Only diagonal components of the dielectric tensor are presented ( xx , yy and zz ); off-diagonals components are close to zero. The Cartesian directions correspond to the coordinate system in Fig. 2.

20

  Fig. 2(a) Fig. 2(b) Fig. 2(c)

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5127513 code

  20 following the method outlined in Ref. 21. The code computes spectral weights w n (k)

	This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.	for d-MAPbI 3 with HEA (without Cs and including relaxation of atomic positions for the entire structure) are presented in SI (Figs. S1 and S2). Overall trends are the same (includ-ing the type-I band alignment, see Fig. 1), but effects are weaker (smaller band gap opening, smaller band discontinuities, and not so large effective mass enhancement). When relaxing positions of Cs atoms, the cutoff energy E cut for the plane-wave expansion was increased by 25% above a recommended maximum value in pseudopotential files (VASP tag PREC=High, E cut = 275 eV). In band structure calculations, E cut corresponded to the recommended maximum value in pseudopotential files (VASP tag PREC=Normal, E cut = 220 eV). The absolute change in the band gap as a result of different E cut is marginal (not more than 5 meV). The Monkhorst and Pack 19 shifted k mesh was used for the structural relaxation (Cs atoms only), while the Γ-centred mesh was selected for calculations of the band gap and optical properties (see SI, Table SI). Linear optical properties were calculated within an independent particle approximation ignoring local field effects. The spin-orbit coupling (SOC) was taken into account. The number of bands was increased twice from its default value when calculating optical properties. Other relevant parameters are summarized in Table SI (see SI). Effective masses were determined at the PBE+SOC level by fitting the energy dispersion E(k) to a second order polynomial within the range of ∆k = ±0.003 Å -1 near to the band edge sampled with 9 -11 intermediate k points. Band structure unfolding and Bloch character calculations were done with the fold2Bloch 3 . Calculations

TABLE I :

 I Band gaps and band alignments (eV) See Fig. 4(a) for definition of band alignments. The positive value of δ v indicates that the VBE is shifted down in energy relative to the stoichiometric phase. The positive value of δ c corresponds to the upward shift of the CBE in energy.

	Parameter		Structure
		stoichiometric	I/Pb = 3.22, x = 0.1	I/Pb = 3.5, x = 0.2
		Fig. 2(a)	Fig. 2(b)	Fig. 2(c)
	E g PBE (no SOC)	1.47	1.59	1.83
	E g PBE (with SOC)	0.38	0.66	1.12
	δ v PBE (with SOC) a	n/a	0.13	0.13
	δ c PBE (with SOC) a	n/a	0.15	0.61

a

TABLE II :

 II Effective masses of electron and holes m * /m e .

This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.PLEASE CITE THIS ARTICLE ASDOI:10.1063/1.5127513

ACKNOWLEDGMENTS

Authors are indebted to Prof. Jacky Even from INSA Rennes for fruitful discussions. C.Z. and O.R. acknowledge funding provided by the Natural Sciences and Engineering Research Council (NSERC) of Canada under the Discovery Grant Program RGPIN-2015-04518. O.R. would like to acknowledge a mobility grant "Materials for green energy" provided by the McMaster University in the framework of International Micro-Fund Initiatives. C.K. and M.K. acknowledge support from Agence Nationale pour la Recherche (MORELESS project). Calculations were performed using a Compute Canada infrastructure supported by the Canada Foundation for Innovation under the John R. Evans Leaders Fund program. C.Z. acknowledges the mobility grant provided by Metropole Rennes. This paper is dedicated to Dr. Jean-François Halet at the occasion of his 60th birthday. It is also a tribute to the memory of Prof. George Papavassiliou.

SUPPLEMENTARY MATERIAL

See the supplementary material for parameters used in electronic structure calculations; structure and band dispersion results for d-MAPbI 3 ; density of states, partial density of states, and the inverse participation ratio analysis. References cited therein are 38,41,42 .


Marie Dalémat

Michel Coret

Adrien Leygue

Erwan Verron

Measuring stress field without constitutive equation

published or not. The documents may come

Introduction

The formulation of deformable mechanical problems involve two types of governing equations. On the one side, the equilibrium (or motion) equations and compatibility conditions (displacement-strain relationships) are uncertainty free; on the other side the constitutive equations which reflect the material response are to use with caution. First, the choice of a constitutive model is rendered difficult by the large number of available models (see for example [START_REF] Marckmann | Comparison of Hyperelastic Models for Rubber-Like Materials[END_REF] for hyperelasticity). Second, the model parameters must be determined thanks to experimental data, which are by nature restricted to specific deformation states. Practically, a constitutive equation can be seen as the way (i) to produce an unnoisy continuous response of a given material from simple mechanical tests, (ii) to interpolate the response between measured data points, but more importantly (iii) to extrapolate to other deformation states, larger strain and/or dimensions (from 1D to 2D or 3D). Obviously such extrapolation is difficult to validate experimentally and may lead to inaccurate numerical results.

The continuous improvement of full-field measurement techniques such as Digital Image Correlation (DIC), presented for example in [START_REF] Sutton | Image Correlation for Shape, Motion and Deformation Measurements: Basic Concepts,Theory and Applications[END_REF], allows the identification of model parameters considering non-standard tests, e.g. experiments with complex geometry and/or loading conditions that involves inhomogeneous fields. Among others, we can mention the Virtual Field Method, the Equilibrium Gap Method, the Constitutive Equation Gap Method, the Reciprocity Gap Method and Finite Element Model Updating; for overviews on theses approaches, the interested reader can refer to [START_REF] Avril | Overview of Identification Methods of Mechanical Parameters Based on Fullfield Measurements[END_REF] and [START_REF] Grédiac | The Virtual Fields Method for Extracting Constitutive Parameters From Full-Field Measurements: a Review[END_REF].

Plus, thanks to the development of "Data Driven" techniques (also referred to as Machine Learning or data mining in various fields, e.g. science, marketing... [START_REF] Berry | Data Mining Techniques: For Marketing, Sales, and Customer Support[END_REF]), new techniques are proposed for materials response in order to assist engineers in both identification and simulation processes. It is to note that unlike Data mining techniques in other fields, physical considerations on material response, e.g. isotropy or other invariances, can greatly help the dedicated solvers [START_REF] Ling | Machine learning strategies for systems with invariance properties[END_REF]. Neural networks are used to determine the material parameters of hyperelastic models [START_REF] Shen | Finite element analysis of V-ribbed belts using neural network based hyperelastic material model[END_REF] and of viscoplastic ones [START_REF] Furukawa | Implicit constitutive modelling for viscoplasticity using neural networks[END_REF]; but also to improve finite element procedures [START_REF] Hashash | Numerical implementation of a neural network based material model in finite element analysis[END_REF]. Manifold learning approaches are considered as well: for example to identify material response using the imprint shapes of indentation tests [START_REF] Meng | Hernot, Identification of material properties using indentation test and shape manifold learning approach[END_REF] and to improve elastodynamic simulations [START_REF] Millán | Nonlinear manifold learning for model reduction in finite elastodynamics[END_REF]. Eventually, the robust identification of the model parameters remains a complex task. Finally, it is sometimes argued that the constitutive equation paradigm has reached a limit; this is the reason why some authors propose to overcome the difficult choice of a constitutive equation by replacing it by the identification of the material response without any underlying constitutive equations. In this different way of thinking, the cornerstone is to develop strainstress databases filled with sufficiently "rich" and "smart" information for simulation tools. Thus, material response must be gathered for a large set of different loading conditions, for which stress is usually difficult to measure. For example, for small strain elastoplasticity, Réthoré et al. measure stress and strain fields from equilibrium equations and DIC measurements by introducing a specific decomposition of the measured strain field and computing the stress fields [START_REF] Réthoré | Computational measurements of stress fields from digital images[END_REF]. Also, Seghir and Pierron demonstrate that stress fields can be computed considering DIC measurements during dynamic experiments [START_REF] Seghir | A Novel Image-based Ultrasonic Test to Map Material Mechanical Properties at High Strain-rates[END_REF]: for samples submitted to ultrasonic vibrations, they measure both displacement and acceleration fields with an ultrahigh speed camera; then the stress field is calculated by considering the motion equation without any force measurement neither constitutive equation. Finally, for hyperelasticity, Latorre et al. propose the concept of What-You-Prescribedis-What-You-Get (WYPiWYG) that consists in replacing the standard strain energy density function with model-free splines [START_REF] Latorre | Experimental data reduction for hyperelasticity[END_REF]. This technique is used with real data; nevertheless it is produced by basic tests where stress fields are easily measured.

In the present paper, a novel method of identification is considered: the "Data-Driven Identification" (DDI) algorithm [START_REF] Leygue | Data-based derivation of material response[END_REF]. It is based on the recent "Data Driven Computational Mechanics" (DDCM) paradigm: bypassing the constitutive equations in mechanical problem formulation [START_REF] Kirchdoerfer | Data-driven computational mechanics[END_REF]. Its originality lies in the absence of an explicit constitutive equation in the governing equations: during finite element computations, the classical constitutive equation is replaced by a large stress-strain database. This theory was first introduced by Ortiz et al. [START_REF] Kirchdoerfer | Data-driven computational mechanics[END_REF] and completed for noise dealing [START_REF] Kirchdoerfer | Data Driven Computing with noisy material data sets[END_REF], large strain [START_REF] Nguyen | A data-driven approach to nonlinear elasticity[END_REF] and dynamic problems [START_REF] Kirchdoerfer | Data-driven computing in dynamics[END_REF]. Roughly speaking, in the DDCM framework the mechanical response of a structure is directly computed considering a large database of "material states", these material states being strain-stress couples that sample the material response. The corresponding algorithm consists in the minimization of the distance between the "mechanical states", i.e. the strain-stress couples to calculate (in each finite element of the mesh for instance) and the material states database, under the constraints of satisfying both equilibrium equations and compatibility conditions. Slightly different formulations of this solver have been proposed [START_REF] Ayensa-Jiménez | A new reliability-based data-driven approach for noisy experimental data with physical constraints[END_REF][START_REF] Kanno | Mixed-Integer Programming Formulation of a Data-Driven Solver in Computational Elasticity[END_REF].

The Data Driven Identification method modifies the DDCM technique for identification of stress fields: for complex experiments, the heterogeneous stress fields are computed without constitutive equation. The success of the method is closely related to the availability of a large database of mechanical problems, i.e. sample geometry and loading conditions and their corresponding displacement fields. This approach has been validated only with synthetic data for truss structures and for small-strain (non linear) elasticity [START_REF] Leygue | Data-based derivation of material response[END_REF]. Here, we aim at validating the DDI approach with real experimental data. Practically, it consists in computing stress fields from displacement fields measured by DIC during the uniaxial extension of a perforated elastomer sheet. In the following, the DDI algorithm is first briefly recalled. Then, the experimental procedure is described; the emphasis is laid on the practical challenges induced by real data as compared to synthetic data: the mechanical problem has to be reformulated at the boundaries due to the unavailability of a few data (close to the grips or holes and in noise affected areas). Once these difficulties overcome, results consists in "measured" heterogeneous stress fields on the membrane surface. The relevance of these results is validated by comparison with standard uniaxial tensile experiments and the experimental derivation of the strain energy density of the material is discussed.

Brief recall of the Data-Driven Identification method

The DDI method consists in computing heterogeneous stress fields without constitutive equation, using a large database of measured displacement fields. In [START_REF] Leygue | Data-based derivation of material response[END_REF], the derivation is proposed for small strain (non linear) elasticity. In the following, only the basics of the method are recalled.

Inputs and outputs

We consider a set of mechanical problems defined by its geometry through a meshed 2D-structure and a given set of loading conditions (denoted by the superscript • Y ). Once the corresponding experiment being conducted, the following quantities are assumed to be available:

• (I-1) the nodal displacements u Y j , j being the node number,

• (I-2) the matrix B Y e j that encodes both geometry and connectivity, e being the quadrature point number. Then, the strain field at point e is defined by (for small strain):

ε Y e = ∑ j B Y e j • u Y j , (1)
• (I-3) the nodal forces f Y j , which are null in the bulk (interior of the mesh) and known on the boundaries;

These quantities are the inputs of the DDI algorithm. The following two are the intrinsic parameters of the methods:

• (In-1) the size N * of the (stress-strain) database that samples the material response,

• (In-2) the norm C that defines the distance between two points in the stressstrain space.

After convergence, the DDI algorithm computes the following outputs:

• (O-1)
for each problem Y , the stress field σ Y e that satisfies the equilibrium equations at each integration point e:

∑ e w Y e B Y e j • σ Y e = f Y j ∀Y, j, (2)
where w Y e is the integration weight of point e. The pair composed by the strain and stress at a given point of the mesh (ε Y e , σ Y e ) is referred to as a mechanical state because it is mechanically admissible, i.e. it fulfills both compatibility and equilibrium equations,

• (O-2) the database of material states (ε * i , σ * i ) that maps the mechanical response of the material. The total number of material states is N * ; the distance to mechanical states is evaluated thanks to the " norm C", which is chosen in an energetic manner, C being a fourth-order positive definite pseudostiffness tensor:

||(ε, σ )|| 2 C = 1 2 (ε : C : ε + σ : C -1 : σ ). (3)

Some details on the DDI method

For a given number of mechanical problems, the algorithm aims at finding the material states that are the closest to the mechanical states, the latest being half known (strain) and half unknown (stress) but constrained by equilibrium equations. Mathematically, it reduces to the following constrained minimization problem:

solution = arg min σ Y e ,ε * e Y ,σ * e Y E (σ Y e , ε * e Y , σ * e Y ) (4)
where

E (σ Y e , ε * e Y , σ * e Y ) = ∑ Y ∑ e w Y e ||(ε Y e -ε * e Y , σ Y e -σ * e Y )|| 2 C , (5)
under the constraints of respecting equilibrium equations [START_REF] Sutton | Image Correlation for Shape, Motion and Deformation Measurements: Basic Concepts,Theory and Applications[END_REF], that the material state (ε * e Y , σ * e Y ) associated to the element e Y belong to the database

(ε * i , σ * i ) N * i=1 .
Practically,

• The initial step computes the material states by the k-means algorithm [START_REF] Macqueen | Some methods for classification and analysis of multivariate observations[END_REF] over the measured mechanical strain and the stresses are set to zero.

• Then, the procedure iteratively 1. computes the material stress to "recenter" them in the cloud formed by the mechanical states (according to the norm C),

2. updates the mechanical stresses to be at equilibrium and as close as possible to the material stress, 3. computes the mapping between mechanical and materials states, 4. updates the material strain, until convergence.

In [START_REF] Leygue | Data-based derivation of material response[END_REF], the validity of the DDI scheme has been assessed thanks to synthetic data.

Application of the DDI method to real data

In the present work, the relevance of the DDI algorithm is demonstrated with real experimental data. For this first "real-life" application, we consider the general context of large strain quasi-static elasticity thanks to the deformation of planar elastomer membranes. First, the experimental procedure is presented; then it is shown that the extension of the previous algorithm to large strain is straightforward. Finally, the main challenges of the application to real data are detailed: they concern the practical difficulties encountered, especially the incompleteness of measurements.

Experimental procedure

Material and samples The elastomer is an unfilled silicone elastomer Bluesil TM RTV 141; it has been chosen because it does not exhibit any Mullins effect and then it can be considered purely elastic [START_REF] Meunier | Mechanical experimental characterisation and numerical modelling of an unfilled silicone rubber[END_REF]. The strain at break, measured with standard uniaxial tensile test, is about 100%. Moreover, the material is assumed homogeneous, isotropic and incompressible. Samples consist in thin rectangular sheets with bulges to fit the grips of the tensile machine; they are cured in an oven accordingly to the technical prescriptions. Once vulcanized, the membrane is perforated with hollow punches. The geometry of the sample is shown in Figure 1. Curing being performed in an open mould, sheets thickness is not prescribed but simply measured afterwards: 3.7 mm ± 0.7mm.

The plane stress assumption is adopted. Coupling it with the incompressibility hypothesis allows the reduction of measurements to planar components: if one Experimental apparatus: loading conditions and measurements The experimental set-up is presented in Figure 2. The uniaxial tensile experiments are conducted with an Instron 3369 electromechanical tensile machine. A specific grip system holds the membrane by clamping the bulges. The loading speed is set to 12 mm/min which corresponds to an overall strain rate of about 1.5•10 -3 s -1 such that experiments can be considered as quasi-static.

The loading force is measured with a force cell. Kinematics are measured thanks to the DIC technique:

• a speckle pattern is applied (by spraying) on the surface of samples (that gives an approximate feature size of 10px),

• a 29 Mpx camera with a 100 mm-focal length lens (the camera is a AlliedVision GT 6600, with an image resolution of 6576 × 4384; the lens is a Zeiss makro-planar T* 2/100), is positionned at a stand-off distance Z ≈ 1 m,

• and an appropriate lighting system is used.

Images are recorded and the displacement field is computed with the commercial software VIC-2D TM (Correlated Solutions). The estimated out-of-plane motion (due to sample thinning) is about ∆Z = 0.3 mm in average over the sample at the last increment, which gives an in-plane strain error of ∆Z/Z ≈ 0.03% according to [START_REF] Sutton | The effect of out-of-plane motion on 2d and 3d digital image correlation measurements[END_REF]. The use of 2D-DIC rather than stereo-DIC is justified. Pictures are recorded every 1.5 s which corresponds to a reasonable engineering strain increment between 2 pictures. The acquisition time of the camera is chosen accordingly to satisfy both correct lighting and no significant change of grey level for one pixel (due to sample motion) during the acquisition. The equivalent pixel size is about 65 µm. The DIC experiments are carried out following 175 the guidelines proposed in [START_REF] Jones | A Good Practices Guide for Digital Image Correlation[END_REF]. An example of recorded a displacement field is presented in Figure 3.

Concerning the order of magnitude of strain, the last loading increment corresponds to about 35% in terms of average engineering strain globally and 100% in terms of local engineering strain. Finally, in order to relate the measurements to the 180 theoretical derivation of Section 2, the different mechanical problems denoted by the superscript Y correspond to the deformation of the samples subjected to different stretching increments. For each of them, the displacement field u Y j is measured with respect to the undeformed configuration.

Extension of the DDI algorithm to large strain

In order to adapt the DDI method to elastomer, Eqs (1)-( 3) must be extended to large strain. For the sake of simplicity, the continuous formulation is chosen in this section, even though the discretised (Finite Element) form of the equations is used in the rest of the paper. Concerning the stress, we adopt an Eulerian description of the problem such that the Cauchy (true) stress σ can be retained. Nevertheless, the corresponding strain tensor must be adapted, the linearized strain tensor ε being not relevant to large strain. Following [START_REF] Holzapfel | Nonlinear Solid Mechanics[END_REF], we adopt the Hencky (true) strain tensor as the natural extension of ε to Eulerian large strain formulation. Considering the deformation gradient F and the left Cauchy-Green strain tensor b = FF T , the Hencky strain tensor is defined by:

ln v = 1 2 ln b, ( 6
)
where v is the (left) pure stretch tensor issued from the polar decomposition of the deformation gradient F = vR. Then, the DDI algorithm is straightforwardly extended to large strain by replacing the linearized strain tensor by the Hencky strain tensor in all the previous equations. As an example, the norm in the stressstrain space Eq. ( 3) becomes:

||(ln v, σ )|| 2 C = 1 2 (ln v : C : ln v + σ : C -1 : σ ). (7)
The original method uses (ε, σ ) and here is introduced the large strain extension. To deal with path-dependency (e.g. plasticity) problems, it might be possible to choose a richer phase space as in [START_REF] Eggersmann | Model-Free Data-Driven inelasticity[END_REF] or [START_REF] Leygue | Nonparametric material state field extraction from full field measurements[END_REF].

Practical challenges when applying the DDI method to "real" data

Application of DDI to real data faces two challenges that are closely related with the discrepancy between the real problem and the theoretical one. It is not possible to obtain mechanical states (strain and stress) at each element because of (i) the incompleteness of measured kinematic fields and (ii) the unavailability of nodal forces but of a unique net force measured through the load cell. In the following, we propose methods to handle these challenges by completing (wisely) the missing mechanical knowledge. These methods are not unique, but our proposals are easily implementable and they are validated on synthetic data in [START_REF] Dalemat | Reliability of the Data-Driven Identification algorithm with respect to incomplete input data[END_REF].

First, the theoretical boundary conditions of the problem must be recalled; they are sketched in Figure 4a. Let us denote Ω the bulk of matter where equilibrium Eq. ( 2) applies simply with f j = 0 (no volume force). The different boundaries are:

• the top boundary Γ T where the load cell gives information. In a perfect (synthetic) case, the equilibrium equation is applied to each node (with nodal forces f j ), • the bottom boundary Γ B where no information is gathered, except displacements,

• the lateral sides of the membrane Γ S . In the perfect case, they correspond to free edges where the mechanical equilibrium is easily assessed,

• the hole edges Γ H for which free edge conditions apply.

Challenge (i): incompleteness of measured kinematic fields

Reasons of incompleteness The incomplete nature of the measured kinematical fields is a classical issue when using DIC [START_REF] Hild | Comparison of Local and Global Approaches to Digital Image Correlation[END_REF]. This full-field measurement technique is based on a local optimization method that consists in computing the displacement vector u which maximizes the correlation between an initial f and a deformed picture g (in grey level), in every subsets defined by few pixels. These subsets are forming the so-called Region Of Interest (ROI) and the correlation is evaluated at each pixel p to find u, such that:

u = arg min ∑ p∈ROI [ f (x p ) -g(x p + u(x p ))] 2 , (8)
where x p is the position of the pixel p in the initial picture. Numerically, u can be approximated thanks to shape functions as in the Finite Element Method (for example see [START_REF] Réthoré | Robust identification of elasto-plastic constitutive law parameters from digital images using 3d kinematics[END_REF]) but non-parametrical approaches are also developed as in VIC-2D TM (the software considered here). In this software, each subset (of N 2 b pixels) furnishes one displacement vector; but an overlap of the subsets is possible (N s is the number of pixels between two subsets). In our case, it leads to one displacement vector u every N s = 7 pixels, based on the computation in a N 2 b = 27 2 -pixel subset (chosen according to both noise and resolution). This naturally forms a structured meshed in which nodal displacements are known every 7 pixels. Using a standard Finite Element approximation, the matrix B Y e j that encodes both geometry and connectivity for the mechanical problem Y , i.e. one deformed image, can be easily computed.

In VIC-2D TM , the accuracy of measurements is controlled by choosing some threshold values above which data are removed. In particular, the "consistency error" checks the back-prediction between neighbour images. It consists in recomputing the supposed initial image from one deformed state image with the inverse displacement fields computed and in checking if the reconstructed initial image and the true initial image are close or not. Moreover, due to measurement noise and/or to large displacement between two successive images, there are some areas in which no displacement is given by the software. Therefore, some data are missing from the measured displacement field.

Definition of new boundaries

When comparing with the theoritical problem shown in Figure 4a, it can be noticed that new boundaries have to be defined for the real problem. They are shown in Fig. 4b : • First, due to inaccurate edge selection during the definition of the ROI (that is manually selected and constrained through N s ), there are some missing zones (and thus nodes) in the neighbourhood of the "true" boundaries: near free surfaces such as holes and sample borders and close to the grips. Thus, it leads to the definition of new boundaries that are close to the theoretical ones and denoted ΓH for hole edges, ΓS for the sides and ΓT for top and ΓB for bottom clamped lines.

• Second, after the DIC process, some data are removed if the consistency error is too large. The corresponding areas are small if displacements of surrounding nodes are well-measured (isolated missing data) or they can spread over few pixels (cluster of missing data). In the first case, the missing displacement value is interpolated from surrounding nodes. In the latter case, the mechanical states cannot be defined because the strain cannot be measured. This leads to the creation of new boundaries denoted Γ C for each cluster of missing data.

Challenge (ii): unavailability of all nodal forces

Practically, it is not possible to access nodal forces that apply on the top clamped line Γ T (or ΓT ): only the total net force F Y cell is measured by the load cell. Therefore, evaluating directly Eq. ( 2) is only possible in the bulk and at free edges where the force is respectively null and null along the normal direction. In our problem, the true free edges are not exactly known as explained in Section 3.3.1. Nevertheless, the side boundaries ΓS can be reasonnably considered as free edges. This is not the case for the hole edges, for which this hypothesis would imply a large error; it has been proven with synthetic data in [START_REF] Dalemat | Reliability of the Data-Driven Identification algorithm with respect to incomplete input data[END_REF] and not recalled here. A practical solution is proposed in the next section.

Completion with mechanical considerations

To summarize, the real problem has been deteriorated compared to the theoretical one, due to the incompleteness of measured kinematic fields (inaccurate edge selection and clusters of missing data) combined with the unavailability of a complete force information. This implies some missing strain or stress and thus impacts the whole mechanical state concerned.

Until then, it has been shown that getting the full mechanical state of the bulk Ω and of the sides ΓS is possible through the measured strain and the computed stress with:

∑ e w Y e B Y e j • σ Y e = 0 for j ∈ ΓFE ∪ Ω. ( 9
)
Note that the free edges assumption is a particular case of the general equilibrium equation: equation σ • n = 0, with the normal vector n, is deductible from the specific geometry of an edge. Therefore, some completion of the mechanical state is proposed thanks to some mechanical considerations:

• The load cell information has no contribution in the algorithm yet. Thus, our belief is that the only way to properly introduce it in the equilibrium equation is to consider that the net force is the sum of non-zero nodal forces f Y j which exerts on ΓT in the loading direction n cell , or so:

∑ j∈ ΓT f Y j • n cell = F Y cell . (10)
Thus, the balance equation on the top boundary is replaced by the force balance in the n cell direction:

∑ j∈ ΓT ∑ e w Y e B Y e j • σ Y e • n cell = F Y cell ∀Y. (11)
Note that we directly consider ΓT instead of Γ T . This approach is less constraining than evenly distributing the net force on the boundary.

• A cluster of missing data, denoted Ω C , cannot be considered in the mechanical problem because no information has been gathered with DIC. We consider the surrounding nodes in which the displacements hava been measured and we denote Γ C the boundary they form. The domain Ω C is not considered in the problem but should be: in Ω C , the matter is mechanically balanced. The equilibrium equations over Ω C Γ C are written and summed. Their collapse in Ω C gives a residual expression of the equilibrium on Γ C :

∑ j∈Γ C ∑ e w Y e B Y e j • σ Y e = 0. (12)
Figure 5: Practical challenges when using the DDI algorithm and our proposals to add mechanical considerations.

This condition represents the overall mechanical balance of the cluster. It can be seen as a zero net force over the boundary Γ C .

• Edges close to the hole cannot be considered as free edges as argued in Section 3.3.2. Yet, the missing matter should be balanced and this case can be treated as previously: a zero net force over the boundary ΓH can be written, similarly than in Eq. ( 12).

In the case of ΓB , there is no need to evaluate the stress at these nodes. Indeed, the global balance brings the stress response with the neighbour nodes in the bulk. A similar way to deal with the boundary consists in treating it as ΓT . Finally, the minimisation process consists in finding the stress closest to the material states under the constraints of mechanical equilibrium (as described in 2.2: step 2 of the algorithm). The constraints are enforced in the algorithm at each nodes for the bulk (Eq. ( 9)) and as a net force over all nodes of the concerned boundary Γ: ∑ j∈Γ f j = F Γ (with f j the mechanical stress contribution at the node j, i.e. the nodal force). The force input F over Γ is the net force over the top boundary (in the sensor direction) (Eq. ( 11)), and is null over the cluster or holes boundaries (Eq. ( 12)).

Graphical summary of the practical challenges

Results and discussion

For the experiments described in Section 3.1, the above-mentionned methods of correction are applied to data issued from DIC; then the DDI algorithm is applied.

Parameters and computation

For the DDI algorithm, the following intrinsic parameters are adopted.

• The number of images, which is here the number of loading increments, is N Y = 165. Thanks to this sufficiently large number of images, the displacement fields measured by DIC do not exhibit too much missing data.

Considering that the DIC process yields meshes with about 8 • 10 5 nodes and N elem ≈ 1.6 • 10 5 linear triangular elements, the total number of computed mechanical states is

N mecha = N elem × N Y ≈ 2.7 • 10 7 .
• The number of material states, which sample the material response, is chosen accordingly to the total number of computed mechanical states:

N * = N mecha 100 ≈ 2.7 • 10 5 .
From a practical point of view, initial values of the material states are computed thanks to the k-means method. The choice of the number of material states N * is also of importance as it changes the mapping between mechanical and material states and so the convergence of the method. A too small number of material states allows the mechanical states (stress in particular) to be far from their associated material state, it adds freedom when solving the algorithm but implies large error. Conversely, a too large number constrains the mechanical stresses to be close to their material data point and the algorithm can stop on local minima. Our previous works on synthetic data suggest one material states every 100 mechanical states.

• The tensor C of the norm Eq. ( 7) is chosen as:

C = C • I (13)
where C = 10 10 Pa and I the fourth-order identity tensor.

Computations are carried out with Matlab on a single machine with 24 cores. Computation time is about 40 h. No particular care has been taken to explicitly parallelize the algorithm. We recall here that the outputs are twofold: all the mechanical states and a database of material states, respectively referred to as (O-1) and (O-2) in Section 2.1. Both spatial distribution and order of magnitude of the stress field are relevant.

Measurement of heterogeneous stress fields

As example, the stress component in the tensile direction σ yy is maximum close to holes highlighting stress concentration and reveals compressive states at top and bottom of these holes. The transversal stress σ xx is noisier: displacements in the x-direction are very small and they lead to insufficiently rich inputs for the DDI algorithm (input (I-1) in Section 2.1) . Moreover, only one mechanical problem is shown; on others, we could observe some white spots in stress distributions: they correspond to clusters of missing data.

Relevance of the results

Discussing the validity of the present approach is not an easy task. Indeed, the naive approach would consist in fitting a hyperelastic constitutive equation, e.g. Mooney-Rivlin or Ogden models, thanks to uniaxial tensile experimental data and comparing stress fields issued from the DDI method with those obtained with Finite Element computations. Nevertheless, such an approach presupposes that numerical results issued from a chosen constitutive model are more valid than the ones measured here. This is not our point of view: we prefer to discuss the results by examining if measured stress fields respect some features of the material response. More precisely, the following questions are investigated: is the isotropy of the material recovered? is its uniaxial tensile response retrieved? In the following, material states (output (O-2)) are used to answer the former question and all measured mechanical states (output (O-1)) to answer the latter one by considering the strain energy density.

Isotropy of the material

The material itself and the moulding process employed for the sample ensure that its mechanical response is isotropic. In the DDI method, material states, i.e. stressstrain couples, are computed with no reference to material symmetry. To verify if measured material states recover isotropy, the misalignment angle between the eigenvectors of strain and stress tensors corresponding to their respective largest principal value is computed for all N * material states. The corresponding results are shown in Figure 7. For more than 90%, respectively 60%, this misalignment is less than ±10.8 • , respectively ±4.6 • . So, using real data, the isotropy of the material is retrieved by the DDI algorithm.

Strain energy density

We consider the uniaxial deformation from the undeformed initial state to a given final mechanical problem Y . This process is parametrized by the time τ such that τ = 0 initially and τ = t for the problem Y . For a given quadrature point defined by its initial position X (in the case of linear triangles, it is the centre of the corresponding finite element) at time τ = 0 and its spatial position x(X, τ) in one of the state in [0,t] the stress power density p int is defined by [START_REF] Holzapfel | Nonlinear Solid Mechanics[END_REF]:

p int (x, τ) = σ (x, τ) : d(x, τ) (14)
where σ (x, τ) is the stress tensor measured by DDI and d(x, τ) is the rate of deformation tensor that can be derived from DIC measurements.

360

Considering now the whole tensile test from 0 to t, the corresponding strain energy density for a given point X is:

w(x(X,t),t) = t 0 p int (x(X, τ), τ) dτ. (15)
By definition, this strain energy is defined per unit of deformed volume. In order to recover the standard strain energy density W (X, τ) defined per unit undeformed volume, we simply invoke the incompressibility of the material such that W (X,t) = w(x(X,t),t).

In the classical isotropic incompressible hyperelasticity theory, the strain energy density is revealed to be a function of the two first strain invariants defined by:

I 1 = tr(b) and I 2 = 1 2 tr(b) 2 -tr(b 2 ) , (17)
the third strain invariant det(b) being equal to 1 due to incompressibility. These two invariants are computed for each point X and time t. strain energy density W (X,t) of all mechanical states, i.e. all elements defined by X for all mechanical problems, obtained using DDI, as a function of these strain invariants.

Finally, the ability to compute the strain energy density can be considered to fit classical hyperelastic models or to use model-free methods such as the What-You-Prescribed-is-What-You-Get (WYPiWYG) method [START_REF] Latorre | Experimental data reduction for hyperelasticity[END_REF].

Remark. In the incompressible framework, the stress field consists in a "material" contribution σ and a "structural" one defined by a hydrostatic pressure p which acts as a Lagrange multiplier of the incompressibility condition:

σ = σ -pI, (18)
I being the 3 × 3 identity tensor. As p depends on the mechanical problem, the stress field σ shown in Fig. 6 does not represent completely the mechanical response of the material. This limitation has been overcome thanks to the strain energy density. Indeed, as the incompressibility constraint can be written as tr(d) = 0, the stress power density Eq. ( 14) reduces to

p int = σ : d. (19)
Thus, the strain energy density issued from Eqs (( 15)-( 16)) describes only the material response.

In order to further discuss our results, we plot the Uniaxial Tensile (UT) strain energy obtained with classical uniaxial tensile tests in Figure 9. At small strain, the accordance is good, but at larger strain, one can observe a mismatch. It is believed that the discrepancy between strain energies from DDI and UT is mainly due to the error on thickness measurement e 0 [START_REF] O'leary | Comparison of methods used to measure the thickness of soft tissues and their influence on the evaluation of tensile stress[END_REF]. For both experiments, the thickness is supposed to be given at ±∆e 0 = ± 0.7mm. This uncertainty is also plotted on the Figure and for graphical clarity both are reported on the uniaxial curve with 2∆e 0 = 1.4 mm.

Conclusion

In this paper, an original coupled experimental-numerical method has been proposed to measure heterogeneous stress fields. Indeed, mixing DIC measurements with the recent Data Driven Identification algorithm allows one to determine the stress field in a complex sample without prescribing any constitutive equation. As a proof of concept, the method has been applied to a multi-perforated elastomer membrane under uniaxial loading conditions. Practical challenges regarding missing data in measured displacement fields and nodal forces have been addressed. DDI algorithm outputs are twofold. First, the mechanical states consisting in heterogeneous stress and strain fields that exactly satisfy equilibrium equations and boundary conditions are derived. Then, the method is revealed able to determine the heterogeneous stress distribution on a planar complex sample using the kinematic DIC measurements. Second, the database of material states that samples the material response is build. Such database is actually the missing piece to the work of Ortiz et al. [START_REF] Kirchdoerfer | Data-driven computational mechanics[END_REF] or Ibañez [START_REF] Ibáñez Pinillo | A manifold learning approach to data-driven computational elasticity and inelasticity[END_REF] to perform mechanical simulation without constitutive equation. Mechanical and material states are relevant, both qualitatively and quantitatively, when compared to classical uniaxial testing data of the material.

Further work is in progress to challenge the method with more complex problems: larger strain, heterogeneous materials, multiaxial loading conditions and inelasticity.

Figure 1 :

 1 Figure 1: Silicone membrane, with bulges and holes. Thickness is about 3.7 mm.

Figure 2 :

 2 Figure 2: Experimental set-up: displacements are gathered with camera and DIC and force measurements with the load cell.

Figure 3 :

 3 Figure 3: Example of a displacement field obtained with DIC at the last increment of the loading path.

Figure 4 :

 4 Figure 4: Fig. a) Theorethical boundaries of the mechanical problem of a perforated membrane under traction. Fig. b) Real boundaries after DIC (incompleteness is due to inaccurate edge selection and clusters of missing data).

Figure 5

 5 Figure 5 summarizes our solutions to the practical challenges.Finally, the minimisation process consists in finding the stress closest to the material states under the constraints of mechanical equilibrium (as described in 2.2: step 2 of the algorithm). The constraints are enforced in the algorithm at each

First, we focus

  on the first output (O-1) of the DDI algorithm: the stress field σ Y e for each mechanical problem. Figure6presents the three in-plane components of the stress field for a given mechanical problem Y .

Figure 6 :

 6 Figure 6: Components of the stress field computed by the DDI algorithm for one mechanical problem. 15

Figure 7 :

 7 Figure 7: Histogram of misalignment angle between principal stress and strain of the material database.

Figure 8

 8 presents the

Figure 8 :

 8 Figure 8: Strain energy density measured by DDI for inhomogeneous experiments with respect to the two first strain invariants.

Figure 9 :

 9 Figure 9: Strain energy densities from the DDI procedure and from uniaxial tensile test, with respect to the two first strain invariants, with uncertainty of thickness measurement.

Acknowledgements This work was performed by using HPC resources of Centrale Nantes Supercomputing Center on the cluster Liger, granted and identified D1705030 by the High Performance Computing Institute (ICI).