Subrata Ghosh

Sujit Kamilya

Titas Pramanik

Mathieu Rouzieres

Radovan Herchel

Sakshi Mehta

Abhishake Mondal

ON/OFF Photoswitching and Thermoinduced Spin Crossover with Cooperative Luminescence in a 2D Iron(II) Coordination Polymer

A 2D coordination polymer, {[Fe(L) 2 (NCSe) 2 ]•6MeOH•14H 2 O} n (1; L = 2,5-dipyridylethynylene-3,4-ethylenedioxythiophene), has been synthesized based on a redox active luminescence ligand. 1 possesses a 2D [4 × 4] square-grid network where the iron(II) center is in a FeN 6 octahedral coordination environment. 1 displays reversible thermoinduced high-spin (HS; S = 2) to diamagnetic low-spin (LS; S = 0) ON/OFF spin-state switching with a T 1/2 value of 150 K. Interestingly, optical reflectivity and photomagnetic studies at 10 K under light irradiation revealed an efficient conversion to a photoinduced metastable HS excited state from a LS ground state. Remarkably, the photoexcited HS state can be reversibly switched ON and OFF by using 625 and 850 nm light-emitting-diode lights. Intriguingly, the thermal dependence of the luminescence intensity of the maximum emission at 524 nm for 1 shows a minimum at around the spin-crossover (SCO) temperature, indicating a cooperative nature between the SCO and luminescence properties. Theoretical calculations confirmed the above findings.

S pin manipulation in molecules triggered by minute alteration of the external stimulus is an attractive and fascinating topic in quantum science and technology. [START_REF] Coronado | Molecular magnetism: from chemical design to spin control in molecules, materials and devices[END_REF][START_REF] Zhang | Magnetic field manipulation of spin current in a single-molecule magnet tunnel junction with two-electron Coulomb interaction[END_REF] Spincrossover (SCO) materials represent a significant way for understanding the molecular spintronics, sensors, switches, nanoscale data recording, and storage devices according to their intrinsic bistable property. [START_REF] Sato | Dynamic molecular crystals with switchable physical properties[END_REF][START_REF] Meng | Manipulating Spin Transition To Achieve Switchable Multifunctions[END_REF][START_REF] Kumar | Emerging trends in spin crossover (SCO) based functional materials and devices[END_REF][START_REF] Kahn | Spin-Transition Polymers: From Molecular Materials Toward Memory Devices[END_REF][START_REF] Halcrow | Spin-crossover materials properties and applications[END_REF] The rational design of a stimuli-responsive multifunctional coordination polymer (CP) by the incorporation of one or more additional cooperative properties, e.g., magnetic, electric, and optical, with the intrinsic properties, e.g., permanent and tunable porosity, versatile structural topologies, and physicochemical properties, has recently gained significant attention. [START_REF] Minguez Espallargas | Magnetic functionalities in MOFs: from the framework to the pore[END_REF][START_REF] Cui | Metal-Organic Frameworks as Platforms for Functional Materials[END_REF][START_REF] Kitagawa | Functional Porous Coordination Polymers[END_REF][START_REF] Halder | Guest-Dependent Spin Crossover in a Nanoporous Molecular Framework Material[END_REF][START_REF] Qiu | Enhanced dielectricity coupled to spin-crossover in a one-dimensional polymer iron(ii) incorporating tetrathiafulvalene[END_REF][START_REF] Rice | Photophysics Modulation in Photoswitchable Metal-Organic Frameworks[END_REF][START_REF] Garg | Conductance Photoswitching of Metal-Organic Frameworks with Embedded Spiropyran[END_REF][START_REF] Feng | Stimuli-responsive multifunctional metal-organic framework nanoparticles for enhanced chemo-photothermal therapy[END_REF][START_REF] Yao | A Dual-Stimuli-Responsive Coordination Network Featuring Reversible Wide-Range Luminescence-Tuning Behavior[END_REF][START_REF] Li | Structure Switching and Modulation of the Magnetic Properties in Diarylethene-Bridged Metallosupramolecular Compounds by Controlled Coordination-Driven Self-Assembly[END_REF] In this aspect, the simultaneous incorporation of both SCO and luminescence properties in a CP is challenging but highly attractive because modulation of the luminescence signal can signify low-spin (LS)/high-spin (HS) population during spin-state switching. [START_REF] Jiao | Synergetic effect between spin crossover and luminescence in the [Fe(bpp) 2 ][BF 4 ] 2 (bpp = 2,6-bis(pyrazol-1-yl)pyridine) complex[END_REF][START_REF] Estrader | A Magneto-optical Molecular Device: Interplay of Spin Crossover, Luminescence, Photomagnetism, and Photochromism[END_REF][START_REF] Benaicha | Interplay between spin-crossover and luminescence in a multifunctional single crystal iron(ii) complex: towards a new generation of molecular sensors[END_REF][START_REF] Yuan | Rhodamine 6G-Labeled Pyridyl Aroylhydrazone Fe(II) Complex Exhibiting Synergetic Spin Crossover and Fluorescence[END_REF][START_REF] Wang | Fluorescence modulation via photoinduced spin crossover switched energy transfer from fluorophores to Fe(II) ions[END_REF][START_REF] Lochenie | Spin-Crossover Iron(II) Coordination Polymer with Fluorescent Properties: Correlation between Emission Properties and Spin State[END_REF][START_REF] Meneses-Sańchez | Extrinsic vs. intrinsic luminescence and their interplay with spin crossover in 3D Hofmann-type coordination polymers[END_REF][START_REF] Delgado | Thermo-and photo-modulation of exciplex fluorescence in a 3D spin crossover Hofmann-type coordination polymer[END_REF] As a successful synthetic strategy to improve cooperativity between SCO and luminescence properties, the direct coupling between the SCO unit and appropriate luminescent organic linkers to form extended structures is significantly important, and it is challenging to interplay such cooperative properties, which have been rarely investigated in 2D multifunctional CPs.

In this aspect, herein we present one 2D CP by introducing Fe(NCSe) 2 with a novel 3,4-ethylenedioxythiophene-based ditopic nitrogen-donor ligand L (Scheme S1), which exhibits cooperative SCO and luminescence properties along with exciting ON/OFF photoswitching behavior.

The reaction of a methanolic solution of Fe:2NCSe with ligand L in dichloromethane (DCM) afforded 1 as red crystals (see Figures S1-S9 for detailed characterization). Variabletemperature X-ray diffraction analyses were performed on a single crystal of 1 at 120 and 200 K, respectively (Tables S1 andS2). 1 crystallizes in the orthorhombic space group Ccca. Its structure consists of the neutral molecular motif [Fe-(L) 2 (NCSe) 2 ] n and interstitial methanol (MeOH) and water solvent molecules (Figures 1 andS10). The overall coordination geometry around each iron center is best described by a slightly distorted octahedron with a FeN 6 environment (CShM program; [START_REF] Casanova | The Rich Stereochemistry of Eight-Vertex Polyhedra: A Continuous Shape Measures Study[END_REF] Table S3), where the iron center is coordinated to the four nitrogen-donor atoms of four monodentate L ligands and two nitrogen-donor atoms of two NCSe -coligands. Ligand L is also acting as a bridging ligand to connect iron centers to form the extended 2D [4 × 4] square-grid networks (Figure 1).

At 200 K, the average Fe-N Py and Fe-N NCSe bond distances are 2.212 and 2.090 Å respectively, which suggested the HS state of the iron center (Figure S11 and Table S2). [START_REF] Sciortino | Two-Dimensional Coordination Polymers with Spin Crossover Functionality[END_REF][START_REF] Mondal | Photoand thermo-induced spin crossover in a cyanide-bridged {Mo V 2 Fe II 2 } rhombus molecule[END_REF][START_REF] Mondal | Photomagnetic effect in a cyanide-bridged mixed-valence {Fe II 2 Fe III 2 } molecular square[END_REF] Upon lowering the temperature to 120 K, the overall structural motif remains unchanged; however, a great difference occurs in the iron(II) coordination sphere. The average Fe-N Py and Fe-N NCSe bond distances reduce to 1.999 and 1.941 Å, respectively (Figure S11 and Table S2), which fall in the range expected for the LS iron(II) ion. [START_REF] Sciortino | Two-Dimensional Coordination Polymers with Spin Crossover Functionality[END_REF][START_REF] Mondal | Photoand thermo-induced spin crossover in a cyanide-bridged {Mo V 2 Fe II 2 } rhombus molecule[END_REF][START_REF] Mondal | Photomagnetic effect in a cyanide-bridged mixed-valence {Fe II 2 Fe III 2 } molecular square[END_REF] The coordinated NCSe -coligands are almost linear [∠N-C-Se = 178.6(7)°]; however, ∠Fe-N-C angles are deviated from linearity [168.4(5)°] at 200 K (Table S2). These values remain almost identical upon lowering of the temperature to 120 K (Table S2). The octahedral distortion parameter Σ value of 12.8 at 200 K compared to the value of 9.9 at 120 K also suggests the HS state of the iron(II) center at a higher temperature (Table S2), [START_REF] Sciortino | Two-Dimensional Coordination Polymers with Spin Crossover Functionality[END_REF][START_REF] Mondal | Photoand thermo-induced spin crossover in a cyanide-bridged {Mo V 2 Fe II 2 } rhombus molecule[END_REF][START_REF] Mondal | Photomagnetic effect in a cyanide-bridged mixed-valence {Fe II 2 Fe III 2 } molecular square[END_REF] supported by CShM analyses (Table S3).

It is worth noting that a contraction in the Fe-N py bond distance of around 0.2 Å is observed, which indicates a complete spin-state switching from the HS to LS state of the iron centers. [START_REF] Sciortino | Two-Dimensional Coordination Polymers with Spin Crossover Functionality[END_REF][START_REF] Mondal | Photoand thermo-induced spin crossover in a cyanide-bridged {Mo V 2 Fe II 2 } rhombus molecule[END_REF][START_REF] Mondal | Photomagnetic effect in a cyanide-bridged mixed-valence {Fe II 2 Fe III 2 } molecular square[END_REF] There is a concomitant decrease in the unit cell volume of around 4.5%, the grid diagonal dimension from {28.55 × 28.22} Å 2 to {28.00 × 27.72} Å 2 and the Fe•••Fe length through the bridging L ligand from 20.07 to 19.70 Å (Figure S12 and Table S2). Around a 6.7% decrease in the solvent-accessible void volume from 4126.6 to 3851.2 Å 3 with a drop from 37.3% to 36.4% of the crystal volume has also been observed.

Magnetic susceptibility measurements were carried out on both fresh crystals kept in a mother liquor (MeOH/DCM/ diethyl ether) and a polycrystalline sample of 1 at a magnetic field of 1000 Oe in the temperature range of 2-280 K (Figure 2). The measured χT [χ is magnetic susceptibility equal to M/ H per iron(II) ion] value at 280 K is 3.25 cm 3 mol -1 K, which lies in the expected range for the iron(II) ion in the HS state (S = 2; g = 2.0; χT = 3.0 cm 3 mol -1 K). [START_REF] Sciortino | Two-Dimensional Coordination Polymers with Spin Crossover Functionality[END_REF][START_REF] Mondal | Photoand thermo-induced spin crossover in a cyanide-bridged {Mo V 2 Fe II 2 } rhombus molecule[END_REF][START_REF] Mondal | Photomagnetic effect in a cyanide-bridged mixed-valence {Fe II 2 Fe III 2 } molecular square[END_REF] Upon lowering of the temperature, the χT value remains almost constant to 204 K and then decreases gradually to reach a value of 0.15 cm 3 mol -1 K at 55 K. The χT value remains almost constant up to 10 K, which indicates a complete and gradual HS-to-LS SCO in 1. [START_REF] Sciortino | Two-Dimensional Coordination Polymers with Spin Crossover Functionality[END_REF][START_REF] Mondal | Photoand thermo-induced spin crossover in a cyanide-bridged {Mo V 2 Fe II 2 } rhombus molecule[END_REF][START_REF] Mondal | Photomagnetic effect in a cyanide-bridged mixed-valence {Fe II 2 Fe III 2 } molecular square[END_REF] In addition, the χT value at 10 K corresponds to less than 2.5% of the residual HS iron(II) center. The small decrease of the χT value below 10 K is probably due to the zero-field splitting of the residual HS iron(II) ion. The field dependence of magnetization studies for 1 were carried out from 0 to 7 T at 1.85, 3, 5, and 8 K (Figure S14). At 1.85 K and 7 T, the magnetization value is 0.11μ B , which is due to the residual HS iron(II) center, matching well with the χT data. In addition, M versus H/T measurement displays nonsuperposition of the magnetic data on a single master curve, probably is due to the presence of magnetic anisotropy coming from the residual HS iron(II) ion. It is worth mentioning that no out-ofphase (χ″) signal under a zero direct-current field of the alternating-current susceptibility was observed in 1, confirming the absence of a slow magnetic relaxation phenomenon (Figure S15).

The peak maxima of the d(χT)/dT versus T plot in the temperature range of 240-80 K indicate the SCO temperature T 1/2 at around 150 K in 1. The identical magnetic responses in both the heating and cooling modes suggest a reversible ON/ OFF SCO for 1. The best fit of the χT versus T data (Figure S16) using the ideal solution model [START_REF] Slichter | Pressure-Induced Electronic Changes in Compounds of Iron[END_REF] gives the thermodynamic parameters ΔH = 27.9(2) kJ mol -1 , T 1/2 = 155(1) K, and ΔS = ΔH/T 1/2 = 180 J K -1 mol -1 , which lie in the range expected for the iron(II) SCO systems. [START_REF] Sciortino | Two-Dimensional Coordination Polymers with Spin Crossover Functionality[END_REF][START_REF] Mondal | Photoand thermo-induced spin crossover in a cyanide-bridged {Mo V 2 Fe II 2 } rhombus molecule[END_REF][START_REF] Mondal | Photomagnetic effect in a cyanide-bridged mixed-valence {Fe II 2 Fe III 2 } molecular square[END_REF] Furthermore, optical-reflectivity studies have also been investigated between 270 and 10 K for 1 (Figure S17). Thermal evolution of the absolute reflectivity clearly indicates SCO between the HS and LS states around 150 K, which is very close to the T 1/2 value.

To explore the photosensitivity of 1, the sample was irradiated with white light (WL; power = 0.5 mW cm -2 ) for ca. 100 min at 10 K. The resultant spectrum significantly differs from the spectra obtained in the dark at both 10 and 270 K (Figure S18, left), suggesting the possibility of formation of the HS iron(II) state from the LS state of 1 at 10 K. The evolution of absolute reflectivity versus time exhibits an increase in the absolute-reflectivity value, attaining a saturation after 20 min of WL irradiation (Figure S18, right). This value is matching well with the one obtained at 270 K, which indicates a complete recovery of the HS state and thus suggests a complete photoinduced SCO of 1 from the LS state at 10 K. After irradiation, the thermal stability of the photoinduced metastable HS state was also investigated by heating the sample in the dark (sweep rate = 0.4 K min -1 ; Figure 3, left), which shows a decrease in the absolute reflectivity due to relaxation of the photoinduced HS state to the diamagnetic LS ground state, where a complete relaxation was observed at around 80 K (T LIESST ). After this temperature, the absolutereflectivity values show a feature similar to that observed in the dark. Additional photoexcitation at 10 K with 14 light-emitting diodes (LEDs; 1050-365 nm) showed that all LEDs were able to populate the paramagnetic HS state (Figure S19), in agreement with the efficiency of WL, [START_REF] Aguila | Switchable Fe/Co Prussian blue networks and molecular analogues[END_REF][START_REF] Ababei | Rational Design of a Photomagnetic Chain: Bridging Single-Molecule Magnets with a Spin-Crossover Complex[END_REF][START_REF] Kamilya | Two-Step Thermoinduced Metal-to-Metal Electron Transfer and ON/OFF Photoswitching in a Molecular [Fe 2 Co 2 ] Square Complex[END_REF] where a maximum efficiency was observed for 625 nm. Accordingly, further photoexcitation studies were carried out with 625 nm LED irradiation at 10 K, which led to similar outcomes (Figures S20 andS21). Interestingly, a weak time relaxation of the photoinduced metastable state in the dark at 10 K can be seen mainly in the first 5 min (Figure S22), which becomes very slow afterward. The effect of the irradiation power was also studied for the most efficient LED, i.e., 625 nm (Figure S22), and when the power is increased, the saturation value is slightly higher and is reached relatively faster, which can be ascribed to a better light penetration at the surface of the sample.

Further optical-reflectivity studies with LED lights were carried out at 10 K by cycling photoexcitation at 625 nm (30 min; power = 16 mW cm -2 ) and a second irradiation with other LEDs to explore the possibility of photoreversibility of 1. The most significant changes in the spectra were observed with the 940-780 nm LEDs, as expected for the HS-to-LS photoconversion (Figure S23) for the iron(II) system. Accordingly, further deexcitation studies were carried out with 850 nm LED irradiation at 10 K. During 850 nm LED irradiation, a drop of the absolute-reflectivity data was observed (Figure S24, left). Interestingly, the resulting spectrum significantly differs from the spectrum obtained after 625 nm LED irradiation at 10 K. It is worth noting that the spectrum obtained after 850 nm LED irradiation is almost similar to the spectrum at 10 K in the dark (Figure S24, right), indicating the possibility of photodeexcitation of the metastable HS state to the LS state for 1. After irradiation, an increase of the temperature produced an augment of the reflectivity due to relaxation of the photodeexcitation state. Such an upturn was observed until 80 K. As expected, the data obtained after this temperature exhibit behavior similar to that in the dark. The HS and LS states can thus be reversibly and efficiently ON/OFF switched by using 625 and 850 nm irradiation over several cycles (Figure 3, right) at 10 K. (Photo)magnetic measurements at 10 K under different diode lasers also confirm the above findings (Figure S25).

To find a correlation between the interesting SCO and luminescence properties in 1, variable-temperature fluorescence measurements were carried out from 270 to 70 K. The ligand L shows an emission band centered at around 565 nm upon excitation with 365 nm at room temperature (Figure S26). 1 exhibit a very interesting luminescence property, showing two broad bands centered at 524 and 570 nm upon excitation at 395 nm (Figure S27), which indicates that the luminescence property for 1 is different from that of the ligand albeit with a decreased intensity. As seen in Figure 4 (left), upon a decrease in the temperature, the luminescence intensity of the maximum emission at 524 nm (I 524 ) initially decreased and then reached the lowest value at around 150 K. Upon a further decrease in the temperature, the emission intensity I 524 started to gradually increase up to 70 K. Minimization of I 524 at around the SCO temperature (T 1/2 ) of 150 K suggested modulation of the luminescence properties with a thermal population of the HS and LS states. The results conclude a remarkable synergistic correlation between the SCO and luminescence properties in 1.

A detailed theoretical CASSCF/NEVPT2 calculation 34 revealed the ligand-field strength (10 Dq LS /10 Dq HS ≈ 1.81) and excitation energies for light-induced excited spin-state trapping (LIESST; 1 Γ → 1 Γ*; Δε LS ≈ 18217 cm -1 ) and reverse LIESST ( 5 Γ → 5 Γ*; Δε HS ≈ 11089 cm -1 ; Figure 4), and the role of ligands on the LIESST and reverse LIESST effects in 1 was studied also by time-dependent density functional theory (Figures S29 andS30), confirming the experimental results.

In conclusion, we have successfully prepared a 2D [4 × 4] square-grid iron(II) CP by integrating Fe(NCSe) 2 with a redox-and luminescence-active ligand. The CP displayed reversible thermo-and photoinduced ON/OFF spin-state switching. Intriguingly, the temperature dependence of the luminescence measurements confirmed a synergic correlation between the SCO and luminescence properties in 1. We believe that this work will deliver an important strategy for the design and development of SCO-and luminescence-based magnetooptical materials and open up new ways for the enhancement of smart switchable materials.

■ ASSOCIATED CONTENT

* sı Supporting Information

The Supporting Information is available free of charge at https://pubs.acs.org/doi/10.1021/acs.inorgchem.0c02136.

Experimental sections, powder X-ray diffraction, thermogravimetric analysis, IR and electronic spectra, crystal packing, magnetic and cyclic voltammetry data, and crystallographic tables (PDF)

Figure 1 .

 1 Figure 1. Coordination environment around the iron(II) (left) and 2D network showing a [4 × 4] square grid (right) in 1. Hydrogen atoms and solvent molecules are omitted for clarity (Fe, green; C, gray; N, blue; O, red; S, orange; Se, brown).

Figure 2 .

 2 Figure 2. Temperature dependence of the χT product for 1 at 1000 Oe in the cooling (blue) and heating (red) modes. Inset: dχT/dT versus T identifies T 1/2 = 150 K for 1.

Figure 3 .

 3 Figure 3. Left: Thermal evolution of the absolute-reflectivity signal plotted at 640 nm in the cooling mode, during WL irradiation (10 K), and in the heating mode for 1. Right: Time evolution of the absolute reflectivity under successive irradiations by 625 and 850 nm LEDs for 1 at 10 K.

Figure 4 .

 4 Figure 4. Left: Thermal variation of the luminescence intensity of the maximum emission (λ em = 524 nm) for 1 with a temperature dependence of the χT product. Middle: Splitting of the d orbitals resulting from AILFT analysis of CASSCF/NEVPT2 for 1. Right: Ligand-field terms resulting from CASSCF/NEVPT2 calculations showing also possible excitation energies to the lowest state with the same multiplicity as the ground state.

  Accession CodesCCDC 2011126 and 2011127 contain the supplementary crystallographic data for this paper. These data can be obtained free of charge via www.ccdc.cam.ac.uk/data_request/cif, or by emailing data_request@ccdc.cam.ac.uk, or by contacting The Cambridge Crystallographic Data Centre, 12 Union Road, Cambridge CB2 1EZ, UK; fax: +44 1223 336033.

	■ AUTHOR INFORMATION Corresponding Author
	Abhishake Mondal -Solid State and Structural Chemistry
	Unit, Indian Institute of Science (IISc), Bangalore 560012,
	India; orcid.org/0000-0002-5061-2326; Email: mondal@
	iisc.ac.in

Downloaded via UNIV OF BORDEAUX on December 14, 2020 at 17:07:48 (UTC).See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

https://dx.doi.org/10.1021/acs.inorgchem.0c02136Inorg. Chem. 2020, 59, 13009-13013


Tiana Deplancke

Marc Fivel

Olivier Lame

D strain rate-dependent constitutive model of UHMWPE: from crystalline network to fibrillar structure behavior

Keywords: constitutive behavior, semicrystalline polymeric material, mechanical testing, microstructures, UHMWPE

    1D strain rate-dependent constitutive model of UHMWPE: From crystalline network to fibrillar structure behavior Tiana Deplancke, Fivel Marc, Olivier Lame

INTRODUCTION

The mechanical behavior of semicrystalline polymers is highly non-linear and sensitive to several external parameters (such as temperature, strain rate) and structural factors (crystallinity, molecular weight, physical entanglements and cross-linking). Semicrystalline polymers may be considered as nano-composite materials. Their structure comprises both amorphous and crystalline phases interconnected in a complex manner.

Over the past two decades several constitutive models were proposed to capture the hyperelasticviscoplastic deformation behavior of solid polymers. Because of the complexity of the problem and due to the lack of a detailed understanding of the local deformation mechanisms, several authors developed phenomenological models [START_REF] Bardenhagen | Three-dimensional, finite deformation, viscoplastic constitutive models for polymeric materials[END_REF][START_REF] Epee | Constitutive model for a semi-crystalline polymer under dynamic loading[END_REF][START_REF] Khan | Finite deformation of a polymer: Experiments and modeling[END_REF][START_REF] Krairi | A thermodynamically-based constitutive model for thermoplastic polymers coupling viscoelasticity, viscoplasticity and ductile damage[END_REF][START_REF] Pouriayevali | A constitutive description of the rate-sensitive response of semi-crystalline polymers[END_REF][START_REF] Yakimets | Model to predict the viscoelastic response of a semicrystalline polymer under complex cyclic mechanical loading and unloading conditions[END_REF].

In the case of amorphous polymers, both glassy and rubbery states have been taken into account in physically based models that gave rather satisfactory and predictive results. They are often inspired by the early work of Haward and Thackray [START_REF] Boyce | Constitutive Models of Rubber Elasticity: A Review[END_REF][START_REF] Buckley | Glass-rubber constitutive model for amorphous polymers near the glass transition[END_REF][START_REF] Haward | The Use of a Mathematical Model to Describe Isothermal Stress-Strain Curves in Glassy Thermoplastics[END_REF][START_REF] Xue | Welding Behavior of Semicrystalline Polymers . 1 . The Effect of Nonequilibrium Chain Conformations on Autoadhesion of UHMWPE[END_REF]. More recently, mechanical modeling connecting both yielding and strain hardening (due to the evolution of the macromolecular network structure) have been also proposed [START_REF] Haward | Strain-hardening thermoplastic[END_REF][START_REF] Senden | Strain Hardening and its Relation to Bauschinger Effects in Oriented Polymers[END_REF]. When considering the semicrystalline polymer both phenomenological [START_REF] Ayoub | Modelling large deformation behaviour under loading -unloading of semicrystalline polymers : Application to a high density polyethylene[END_REF][START_REF] Drozdov | Mullins' effect in semicrystalline polymers: Experiments and modeling[END_REF][START_REF] Dusunceli | Modelling effects of degree of crystallinity on mechanical behavior of semicrystalline polymers[END_REF][START_REF] Khan | Thermo-mechanical large deformation response and constitutive modeling of viscoelastic polymers over a wide range of strain rates and temperatures[END_REF][START_REF] Lee | A thermo-mechanical large deformation constitutive model for polymers based on material network description: Application to a semi-crystalline polyamide 66[END_REF][START_REF] Regrain | Multi-mechanism models for semicrystalline polymer: Constitutive relations and finite element implementation[END_REF] and micromacro approaches have been proposed [START_REF] Dommelen | Micromechanical modeling of the elasto-viscoplastic behavior of semi-crystalline polymers[END_REF][START_REF] Gueguen | Mechanics of Materials A new threephase model to estimate the effective elastic properties of semi-crystalline polymers : Application to PET[END_REF][START_REF] Nikolov | Multi-scale constitutive modeling of the small deformations of semi-crystalline polymers[END_REF][START_REF] Sedighiamiri | A micromechanical study on the deformation kinetics of oriented semicrystalline polymers[END_REF][START_REF] Shojaei | Viscoplasticity analysis of semicrystalline polymers : A multiscale approach within micromechanics framework[END_REF][START_REF] Uchida | Finite element simulation of deformation behavior of semi-crystalline polymers with multi-spherulitic mesostructure[END_REF].

All these approaches can be accurate in a certain domain, but their main drawback is that they cannot really capture the huge microstructural evolution that occurs during fibrillation of a semicrystalline polymer. As a consequence, their reliability beyond their fitted domain remains limited.

On the other hand, physical approaches (generally 1D models) that have been developed since the 70s have allowed establishing the main deformation mechanisms that occur during the elastic regime, the fibrillation stage and the fracture process. Viscoelasticity, plasticity, damage, hardening and fibrillation fractures with crack propagation are qualitatively well described in the literature [START_REF] Bouvard | Review of Hierarchical Multiscale Modeling to Describe the Mechanical Behavior of[END_REF][START_REF] Michler | Mechanical Properties of Polymers based on Nanostructure and Morphology[END_REF][START_REF] Ward | The mechanical properties of solid polymers[END_REF].

Since macroscopic approaches do not take into account the huge microstructural evolution occurring during fibrillation (crystalline network collapse), plasticity is generally viewed as a continuous process from yielding to hardening. This is one of the main weak points of semi crystalline modeling that should be addressed in a physical approach as the mechanisms of fibrillation are qualitatively well understood, and can accurately reproduce a memory effect between the initial and fibrillar structure [START_REF] Xiong | On the strain-induced fibrillar microstructure of polyethylene: Influence of chemical structure, initial morphology and draw temperature[END_REF]. Basically, plastic flow leads to strong softening due to fractures of crystalline lamellae into small blocks. However, entanglements and molecular intercrystalline connections progressively catch up this collapse leading finally to hardening. It is largely agreed that these connections are essentially contained in the confined amorphous phase and behaves as Stress Transmitters (ST see [START_REF] Seguela | Critical review of the molecular topology of semicrystalline polymers: The origin and assessment of intercrystalline tie molecules and chain entanglements[END_REF]). They are supposed to be essentially composed of chains connecting adjacent crystallites, namely, Tie Molecules (TMs) which are known to strongly affect the mechanical properties of semicrystalline polymers both in the elastic domain as well as during the plastic behavior. Their direct quantification remains however impossible. According to the work of [START_REF] Huang | Dependence of slow crack growth in polyethylene on butyl branch density: Morphology and theory[END_REF], the probability of forming a tie molecule increases with the molecular mass of the polymer. As a consequence, it exists a critical molecular weight above which chains becomes intertwined into a macromolecular network [START_REF] Huang | Dependence of slow crack growth in polyethylene on butyl branch density: Morphology and theory[END_REF][START_REF] Porter | The entanglement concept in polymer systems[END_REF][START_REF] Seguela | Critical review of the molecular topology of semicrystalline polymers: The origin and assessment of intercrystalline tie molecules and chain entanglements[END_REF]. Some modeling approaches have been proposed to establish connections between the macroscopic mechanical behavior and this specific amorphous phase [START_REF] G'sell | Evolution of microstructure in semi-crystalline polymers under large plastic deformation[END_REF][START_REF] Lee | A thermo-mechanical large deformation constitutive model for polymers based on material network description: Application to a semi-crystalline polyamide 66[END_REF][START_REF] Oshmyan | Simulation of small-strain deformations of semicrystalline polymer: Coupling of structural transformations with stress-strain response[END_REF][START_REF] Roguet | Mechanical features of the rubbery amorphous phase in tension and torsion in a semi-crystalline polymer[END_REF].

In the case of Ultra High Molecular Weight Polyethylene (UHMWPE), the length of the chains implies that nearly all the chains are tie chains, and this asserts a very strong macromolecular network that is generally used to explain the unusual capabilities to accumulate deformation without necking [START_REF] Bonten | A New Hypothesis To Describe the Mechanisms Acting[END_REF][START_REF] Humbert | Polyethylene yielding behaviour: What is behind the correlation between yield stress and crystallinity?[END_REF]. Very few attempts exist for modeling the mechanical behavior of such a polymer [START_REF] Bergström | Constitutive modeling of ultrahigh molecular weight polyethylene under large-deformation and cyclic loading conditions[END_REF], and microstructural specific features are not explicitly taken into account. The strong macromolecular network should clearly be included in a realistic mechanical constitutive law.

UHMWPE is also known to exhibit exceptional wear and impact resistance. As an example, a recent study has highlighted the outstanding shock resistance of UHMWPE in the case of hydrodynamic cavitation erosion (Deplancke et al., 2015a). However, it remains very difficult to predict the plastic deformation and the life time of the coatings under such severe conditions and thus numerical modeling is required. This implies to set up an accurate model of elasto-viscoplastic behavior of UHMWPE capable to take care of large strain rates. This paper presents a mixt approach of the mechanical behavior of UHMWPE combining experimental investigations of sintered UHMWPE with various molecular weights under high strain rate (comparable to impact induced by cavitation erosion) and numerical modeling of the same materials for a large range of strain rates.

To predict the behavior of UHMWPE we propose first a general physical approach aiming at describing the mechanical behavior of the semicrystalline polymers both in the elastic and plastic regime by including the main physical mechanisms existing in the literature. Then, the model will be adapted to UHMWPE by adjusting some of its parameters, in particular to take into account the huge strain recovery observed on these materials even after a large plastic deformation. Such a model could ultimately be used to explain the resistance to cavitation erosion as supposed to happen at high strain rate [START_REF] Roy | Determination of cavitation load spectra -Part 1: Static finite element approach[END_REF].

EXPERIMENTAL PART

Materials

The UHMWPE under investigation is elaborated from three different commercial nascent reactor powders produced by Celanese (ex-Ticona, Oberhaussen, Gernamy). Their average molecular weights deduced from viscosity measurements are ranged between Mv = 0.6, 3.9 and 10.5 million g/mol (GHR 8110, GUR 4113 and 4170 respectively). This nascent powder was processed by sintering performed in an Instron compressive testing machine equipped with an oven "Servantin".

The sintering protocol was optimized in order to obtain the best mechanical properties (see (Deplancke et al., 2015b[START_REF] Deplancke | Diffusion versus cocrystallization of very long polymer chains at interfaces: Experimental study of sintering of UHMWPE nascent powder[END_REF] for details of the protocol).

Cylindrical compression tests specimens with 8 mm in diameter and 8-9 mm in height were then machined out from the sintered samples.

The crystallinity ratio was measured by DSC using a Perkin Elmer Pyris apparatus at a heating rate of 10 °C/min, and calculated considering that the melting enthalpy of a perfect PE crystal is , 2003).

𝐻 𝑚 0 = 290 J/g (Wunderlich

Mechanical characterization

Compression tests have been performed on the three materials with different molecular weights, at strain rates ranging between 10 -4 and 10 s -1 using a conventional servo-hydraulic MTS compression testing apparatus at room temperature with a load cell of 100 kN and compression platens with ball-joint system.

To reduce the friction at the platen/specimen interface, the platens were polished and thin Teflon (PTFE) films were used as lubricating layers between the specimen and the platens.

These compression tests have been complemented by Split Hopkinson Pressure Bar (SHPB) tests for accessing higher strain rates of 10 3 s -1 . This experimental device is based on wave propagation and the specimen mechanical impedance must be adapted to the bars impedances [START_REF] Gama | Hopkinson bar experimental technique: A critical review[END_REF]. In the case of UHMWPE, 2 m long and 38.2 mm diameter magnesium alloy Mg ZK60 T5

bars have been used. The gas gun is filled with compressed air at a pressure of around 4 bars, which represents a mean strain rate of 10 3 s -1 . The three-wave technic [START_REF] Gama | Hopkinson bar experimental technique: A critical review[END_REF]) is used to obtain stress strain curves (see [START_REF] Pfeiffer | Analyse et modélisation du fraisage du bois vert[END_REF] for details of the protocol and device).

Assuming homogeneous deformation of all samples irrespective of drawing conditions, the true strain was computed as ε = ln (L/L0), where L0 is the initial sample gauge length and L its length under straining. The true stress was defined as σ = Fλ /S0, where F is the tensile force, λ = L/L0 is the draw ratio, and S0 is the initial sample cross section. This computation assumes isovolumic strain, neglecting Poisson's effect and cavitation.

RESULTS AND ANALYSIS

The figure 1 displays selected results of compression tests for different strain rates between 10 -4 and 10 3 s -1 on the sintered UHMWPE with 3 different molecular weights 0.6, 3.9 and 10.5 Mg/mol.

Each true stress-strain curve displayed on these plots is a representative curve averaged among four realizations of independent compressive tests.

These curves show:

-A dispersion of the elastic modulus certainly due to the measurement made without any extensometer. To avoid this problem, the data used for numerical identification will be the average elastic moduli of the four compressive tests, measured between 1 and 3% deformation for each strain rate.

-An increase in yield stress with the strain rate followed by strain hardening. The slope of the strain hardening increase also with the strain rate.

-An immediate recovery, of the order of 20%, is observable during the unload stage which was performed at the same strain rate as the loading part. This immediate recovery is strongly dependent of the molecular weight.

Recovery tests over long periods of time, of the order of the month, display a recovery close to 100% both for the molecular weight of 10.5 Mg/mol and 3.9Mg/mol and 85% for the molecular weight of 0.6 Mg/mol in the case of the strain rates between 10 -4 s -1 and 1 s -1 . Indeed, measurements of specimen dimensions after these compression tests reveal almost no change in shape and a length equal to the initial length is recovered in the case of the molecular weight of 10.5 Mg/mol after one month of relaxation.

Figure 1: Experimental strain-stress curves from compressive tests on the three molecular weights and various strain rates in the range from 10 -4 s -1 to 1000s -1 .

Figure 2 shows typical compression tests realized on the sintered UHMWPE for 2 different maximum deformations (40% and 80%) for the molecular weights of 0.6 and 10.5 Mg/mol and at a strain rate of 0.1 s -1 .

The first striking observation is that the shape of the 'unload part' of the curves are very similar regardless of the maximal deformation. Observation of the deformed samples show a barelling effect for plastic deformation above 0.5 so that it could explain the change of slope visible on the curves at this cumulated true strain. This artefact comes from the estimation of the true strain based on a non-local measuremnt of the elongation. We will not consider this change of slope in our modeling.

Consequently, the model presented in this article will take into account the following experimental observations of the mechanical behavior of sintered UHMWPE:

-First a recovery in very short times occurs whatever the speed of deformation and the form of this recovery is very similar for several deformations.

-In addition, one of the specificities of the UHMWPE is that a recovery in longer times (of the order of the month) allows a quasi-complete recovery for deformation as large as 40%.

Therefore, the model proposed in this article should allow to reproduce these two experimental observations and should also remain in agreement with all the observations of the literature concerning the mechanical behavior of the semicrystalline polymers.

CONSTITUTIVE MODEL OF UHMWPE

The main purpose of this part is to propose a uniaxial constitutive model able to reproduce (at least semi-quantitatively) the hyperelasto-viscoplastic response of the three UHMWPE materials for various molecular weights between 0.6 and 10.5 Mg/mol over a wide range of strain-rates: from 10 -4 to 10 3 s -1 . As a uniaxial model, it is relatively easy to include all the physical mechanisms on which there is a consensus in order to i) obtain a model which is satisfactory and understandable from a physical point of view and ii) to promote a better reliability when varying parameters such as the strain rate, crystallinity, etc.

As soon as the crystallinity exceeds a few percent, it is generally admitted that a semicrystalline polymer can be considered as a macromolecular network (a rubber when T>Tg) bridled by a crystalline network formed by lamellar stacking (Detrez, F;Séguéla, R.;[START_REF] Detrez | Quantification et modelisation du couplage plasticite endommagement dans les polymeres semicristallins a morphologie spherolitique. 19 ème congrès Français Mécanique[END_REF][START_REF] Detrez | Plasticity / damage coupling in semi-crystalline polymers prior to yielding : Micromechanisms and damage law identi fi cation[END_REF] generally arranged in an isotropic spherulitic structure see figure 3a. It is then reasonable to consider the decomposition of the global stress into two terms corresponding to the two superimposed networks (see figure 3). In the case of medium molecular mass PE, the macromolecular network consists essentially of the confined amorphous phase which can flow relatively easily at high stress. This assumption is close to the one proposed by Ward and Ayoud [START_REF] Ayoub | Effects of crystal content on the mechanical behaviour of polyethylene under finite strains : Experiments and constitutive modeling[END_REF][START_REF] Ayoub | Modelling large deformation behaviour under loading -unloading of semicrystalline polymers : Application to a high density polyethylene[END_REF][START_REF] Ward | The Role of Molecular Networks and Thermally Activated Processes in the Deformation-Behavior of Polymers[END_REF] and used in several modeling approaches.

In the case of UHMWPE, short range macromolecular network still exists and is even reinforced, but it should also exist a long range macromolecular network that strongly prevents the material from flowing at high stress. Indeed, one should keep in mind that for UHMWPE even after a very high deformation, a nearly full recovery can be observed after several weeks.

As a consequence, we propose first to decompose the stress in two parts: the macromolecular network and the crystalline stack network. The first is strongly involved in the hardening even at low strain and, of course, the long-term strain recovery. It is particularly important for UHMWPE.

The second represents what is common to all semicrystalline polymers i.e. the crystalline stack mechanical behavior. It has to be noticed that this part contains the confined amorphous phase, so that this modeling is itself splitted in two different parts. For large strain, the crystalline network has collapsed and reorganized into small crystalline blocks that are positioned relatively regularly along a fibril. This collapse induces a large decrease of the material stiffness whereas a progressive mechanical hardening is simultaneously observed. The latter is due to the strong alignment of the chains. Indeed, for a large strain, the macromolecular network is strongly stretched which compensates for the loss of mechanical properties of the crystalline network (see figure 4).

Description of the macromolecular network branch

In order to model the stretching of the macromolecular network within the strain-hardening stage, the eight-chain network model, developed by Arruda and Boyce [START_REF] Arruda | A three-dimensional constitutive model for the large stretch behavior of rubber elastic materials[END_REF] is selected for its physical consistency.

𝜎 𝑀𝑁 = 𝐶 0 (𝜆 𝑀𝑁 2 -1 𝜆 𝑀𝑁 ⁄ ) = 𝜇 𝑀𝑁 3 √𝑁 𝑀𝑁 𝜆 𝑐ℎ𝑎𝑖𝑛𝑀𝑁 ℒ -1 ( 𝜆 𝑐ℎ𝑎𝑖𝑛𝑀𝑁 √𝑁 𝑀𝑁 ) (𝜆 𝑀𝑁 2 -1 𝜆 𝑀𝑁 ⁄ ) (2)
In equation ( 2), 𝜆 𝑀𝑁 = 𝜀 𝑀𝑁 + 1 is the draw ratio of the macromolecular network,

𝜆 𝑐ℎ𝑎𝑖𝑛𝑀𝑁 = √𝜆𝑀𝑁 2 -2 𝜆 𝑀𝑁 ⁄ √3
is the stretch of a chain in the eight-chain network in the case of onedimensional mechanical testing, ℒ -1 is the inverse Langevin function, √𝑁 𝑀𝑁 the limiting chain extensibility of the macromolecular network and 𝜇 𝑀 the rubbery modulus.

This model however does not capture the increase of the strain-hardening modulus, with the strain rate.Following the work of Govaert and Haward [START_REF] Goveart | Does the Strain Hardening Modulus of Glassy Polymers Scale with the Flow Stress ?[END_REF][START_REF] Haward | Strain-hardening thermoplastic[END_REF][START_REF] Senden | Strain Hardening and its Relation to Bauschinger Effects in Oriented Polymers[END_REF], a viscous contribution is added to the entropic elasticity. This flow stress is related to the intermolecular frictions at a monomeric or segmental scale. The viscous flow stress is a s,e s,e function of the strain rate 𝜀, which is described by an Eyring relation expressed as a function of the total deformation:

𝐶 𝑎𝑝𝑝 (𝜀) = 𝑘 𝐵 𝑇 𝛼 𝑀𝑁 𝑠𝑖𝑛ℎ -1 ( 𝜀̇ ℎ 𝑘 𝐵 𝑇 exp ( ∆𝐻 𝑀𝑁 𝑘 𝐵 𝑇 )) (3)
With 𝜀̇ the strain rate, ∆𝐻 𝑀 the activation energy and 𝛼 𝑀 the activation volume.

Therefore, the elastic stress in uniaxial loading is represented by the following relation:

𝜎 𝑀𝑁 = (𝐶 0 + 𝐶 𝑎𝑝𝑝 (𝜀)) (𝜆 𝑀𝑁 2 -1 𝜆 𝑀𝑁 ⁄ ) (4)
where 𝐶 0 isthe Boyce modulus, 𝐶 𝑎𝑝𝑝 is the viscous contribution of the strain hardening and 𝜆 𝑀𝑁 = 𝜀 + 1 is the draw ratio of the macromolecular network.

A similar elastic model has been successfully applied to semicrystalline polymers [START_REF] Haward | Strain-hardening thermoplastic[END_REF] and it was able to capture the essential of the physics of the hardening behavior. The model can be schematically represented by two springs, one is strain rate dependent whereas the second is not However, even in the region of very large strain, this phenomenon should be weak for most of the UHMWPE when the temperature remains below the melting point. Thus, the flow behavior of the Version acceptée pour publication dans Mechanics of Materials 137 (2019) 103129 14 macromolecular network is described by a viscosity proportional to the cubic power of the molecular weight as already proposed by [START_REF] Doi | The theory of polymer dynamics[END_REF]Edwards, 1988 andWang and[START_REF] Wang | On the viscosity-temperature behavior of polymer melts[END_REF] :

𝜂 𝑓𝑙𝑜𝑤 = 𝐺 𝑓𝑙𝑜𝑤 𝑀 3 (5)
where M is the molecular weight and 𝜂 𝑓𝑙𝑜𝑤 the viscosity of the macromolecular network.

Consequently, the plastic deformation of the macromolecular network is given by:

𝜀ḟ 𝑙𝑜𝑤 = 𝜎 𝑀𝑁 𝜂 𝑓𝑙𝑜𝑤 (6)
For most of the UHMWPE, M being huge, this viscosity can be considered (in a first approach) as infinite.

Description of the Crystalline network branch

Two springs, relative to crystalline and confined amorphous parts are used to represent the initial modulus of the material. A simple analytical model with an adjustable mechanical coupling is used to take into account the fact that these two phases are neither in parallel nor in series but rather a mixture of both as expressed in equation ( 7):

𝜎 𝐶𝑟𝑁 = 𝜑𝜎 𝑐𝑟 + (1 -𝜑)𝜎 𝑐𝑜𝑛𝑓𝑖𝑛𝑒𝑑 𝑎𝑚𝑜𝑟𝑝ℎ𝑜𝑢𝑠 𝑝ℎ𝑎𝑠𝑒(𝐶𝐴) = (1 -𝜑 + 𝑘𝜑)𝜎 𝐶𝐴 . ( 7
)
where 𝜑 is the crystallinity obtained by DSC experiments and 𝑘 = 𝜎 𝐶𝑟 /𝜎 𝐶𝐴 the adjustable mechanical coupling parameter where 𝜎 𝑐𝑟 and 𝜎 𝐶𝐴 are respectively the stress in the crystal and in the confined amorphous phase. parallel coupling. This model is close to a Halpin Tsai model [START_REF] Halpin | The Halpin-Tsai equations: A review[END_REF], and is used mainly for the sake of simplicity and because it gives reliable results for medium to high density polyethylene including UHMWPE. It is also suitable to take into account the crystal network collapse occurring during plasticity (see the general description of the constitutive model).

The initial value of the parameter k is fitted on the evolution of the Young's modulus as a function of the crystallinity and it evolves during the plastic deformation of the material. The stress on the crystal network is given by 𝜎 𝑐 = 𝐸 𝑐 𝜀 𝐶 with Ec = 28GPa [START_REF] Boyd | The modulus of the amorphous component in polyethylenes[END_REF].

To extrapolate the confined amorphous phase behavior for larger deformation, i.e.during plastic strain, a hyperplastic model has been proposed by [START_REF] Arruda | A three-dimensional constitutive model for the large stretch behavior of rubber elastic materials[END_REF].

𝜎 𝐶𝐴 = 𝜇 𝐶𝐴 3 √𝑁 𝐶𝐴 𝜆 𝑐ℎ𝑎𝑖𝑛𝐶𝐴 ℒ -1 ( 𝜆 𝑐ℎ𝑎𝑖𝑛𝐶𝐴 √𝑁 𝐶𝐴 ) (𝜆 𝐶𝐴 2 -1 𝜆 𝐶𝐴 ⁄ ) (8)
In equation ( 8), 𝜆 𝐶𝐴 = 𝜀 𝐶𝐴 + 1 is the draw ratio,

𝜆 𝑐ℎ𝑎𝑖𝑛𝐶𝐴 = √𝜆𝐶𝐴 2 -2 𝜆 𝐶𝐴 ⁄ √3
is the stretch of a chain in the eight-chain network in the case of one-dimensional mechanical testing, ℒ -1 is the inverse Langevin function, √𝑁 𝐶𝐴 the limiting chain extensibility of the tie molecules which is strongly weaker than that of the rubbery network, and 𝜇 𝐶𝐴 the rubbery modulus.

The amorphous phase confined in the crystallites stacking is not directly comparable to a bulk rubber because of the presence of an interphase, and because they are locked between two crystallites. We have shown that for particular orientations of the stacking, the bulk modulus of the amorphous phase can be more important than the tensile modulus [START_REF] Xiong | Amorphous phase modulus and micro-macro scale relationship in polyethylene via in situ SAXS and WAXS[END_REF]. As a consequence, in some cases the apparent tensile modulus of the amorphous phase could be close to 300 MPa. Consequently, the apparent average modulus of the amorphous phase in a semicrystalline polymer has nothing to do with the one of a pure rubber and is much higher. Using our experimental data, an initial value of 70 MPa seems to be appropriate.

Moreover, Tie Molecules (TM) which link adjacent crystallites are also partially responsible for such a high modulus of the amorphous phase [START_REF] Xiong | Amorphous phase modulus and micro-macro scale relationship in polyethylene via in situ SAXS and WAXS[END_REF]. Indeed, it has been shown that they are initially taut, explaining why they can affect elastic domain and yielding behavior, as well as fibrillar formation [START_REF] Men | Role of the Entangled Amorphous Network in Tensile Deformation of Semicrystalline Polymers[END_REF][START_REF] Seguela | Critical review of the molecular topology of semicrystalline polymers: The origin and assessment of intercrystalline tie molecules and chain entanglements[END_REF][START_REF] Xiong | On the strain-induced fibrillar microstructure of polyethylene: Influence of chemical structure, initial morphology and draw temperature[END_REF] and strain hardening.

Similarly to a rubber, these chains can be considered as active chains [START_REF] Men | Role of the Entangled Amorphous Network in Tensile Deformation of Semicrystalline Polymers[END_REF][START_REF] Patel | Theoretical prediction of tie-chain concentration and its characterization using postyield response[END_REF][START_REF] Peterlin | Folded chain concept of fiber structure[END_REF][START_REF] Seguela | Critical review of the molecular topology of semicrystalline polymers: The origin and assessment of intercrystalline tie molecules and chain entanglements[END_REF]. Due to its entropic elasticity, this network plays a major role in the spring back effect of the material as they still connect the crystalline blocks after fibrillation. Hyper-elastic model is then suitable to describe its behavior considering a very weak limit of extensibility.

Plastic behavior

The plasticity in semicrystalline polymers is known to be due to the motion of dislocations in the crystalline phase which is a thermally activated phenomenon. It is also known that the energy barrier for dislocation motion is increasing with the crystal thickness, and then with the crystallinity.

Based on thermally activated mechanisms and stress driven considerations of defect mobility, the Eyring models [START_REF] Eyring | Viscosity, Plasticity, and Diffusion as Examples of Absolute Reaction Rates[END_REF][START_REF] Gueguen | Micromechanically based formulation of the cooperative model for the yield behavior of semi-crystalline polymers[END_REF][START_REF] Richeton | A formulation of the cooperative model for the yield stress of amorphous polymers for a wide range of strain rates and temperatures[END_REF][START_REF] Ward | The Role of Molecular Networks and Thermally Activated Processes in the Deformation-Behavior of Polymers[END_REF] are classically used to describe the plastic behavior of polymers. For semicrystalline polymers one of the simplest Eyring model can be derived as follows:

𝜀ṗ 𝑙𝑎𝑠𝑡𝑖𝑐 = 𝑘 𝐵 𝑇 ℎ 𝑒 -∆𝐻 𝑃𝑙 𝑘 𝐵 𝑇 ⁄ sinh ( 𝛼 𝑃𝑙 𝜎 𝑘 𝐵 𝑇 ) (9)
with 𝑘 𝐵 the Boltzmann's constant, T the temperature (here equal to the room temperature 293 K), h the Planck constant. Parameters ∆𝐻 𝑃𝑙 and 𝛼 𝑃𝑙 , correspond respectively to the activation energy and the activation volume and could be dependent on the crystallinity. This model explicitly predicts an increase in the yield stress with the increase of the strain rate.

During plastic deformation, following the Peterlin scheme [START_REF] Bartczak | Plasticity of Semicrystalline Polymers[END_REF][START_REF] Peterlin | Plastic deformation of polymers with fibrous structure[END_REF][START_REF] Peterlin | Folded chain concept of fiber structure[END_REF], the microstructure evolves from spherulitic organization where the crystal network is percolating toward a fibrillar microstructure where crystalline blocks are nearly isolated and embedded in the amorphous phase. During this process, the stress concentration in the crystal phase decreases. As a consequence, the mechanical coupling k in eq.7 evolves from the initial situation where it is intermediate between series and parallel toward a situation where the coupling can be assumed to be close to a series coupling i.e. k = 1. The fibrillar structure is supposed to be achieved for a deformation associated to the natural draw ratio [START_REF] Gent | Plastic Yielding of Partially Crystalline Polymers[END_REF][START_REF] Peterlin | Drawing and extrusion of semi-crystalline polymers[END_REF][START_REF] Seguela | On the Natural Draw Ratio of Semi-Crystalline Polymers: Review of the Mechanical, Physical and Molecular Aspects[END_REF][START_REF] Ward | The Role of Molecular Networks and Thermally Activated Processes in the Deformation-Behavior of Polymers[END_REF].

As a consequence, we propose a strong connection between eplastic and k in the crystalline network branch. The following empirical equation is proposed to describe a progressive collapse of the crystalline network as the plastic deformation increases:

𝑘 = 𝑘 𝑖𝑛𝑖𝑡 𝑒 -𝜀 𝑃𝑙 𝜀 𝑁𝑎𝑡_𝐶𝑟 Ln 𝑘 𝑖𝑛𝑖𝑡 (10)
Where 𝑘 𝑖𝑛𝑖𝑡 is determined by fitting a wide range of polyethylene with varying their crystallinity [START_REF] Humbert | A re-examination of the elastic modulus dependence on crystallinity in semi-crystalline polymers[END_REF], 𝜀 𝑃𝑙 is determined from Eyring's expression (eq 3), and 𝜀 𝑁𝑎𝑡_𝐶𝑟 is the deformation for which the fibrillar state is reached (corresponding to the series coupling and consequently k =1).

In other words, when the fibrillary state is attained, the rheological model of the crystal branch transforms into the series model as shown in the following scheme (fig 7). The form of this equation remains phenomenological, but could be justified for low and medium density polyethylene where the first breakages in the crystalline network are not critical, leading finally to a progressively accelerated decrease of parameter k.

At larger deformation, it has been experimentally observed that secondary plasticity could occur in the necked region and could even lead to a second necking formation. This plasticity is associated with the slipping of the fibrils or break of the tie molecule and/or partial disentanglement. It could also be unfolding mechanisms or chains pull out [START_REF] Castagnet | Correlation between macroscopical viscoelastic behaviour and micromechanisms in strained α polyvinylidene fluoride (PVDF)[END_REF] in the crystal blocks.

To capture this plasticity a second non-linear dashpot is added in a series coupling with the spring representative of the confined amorphous phase using the following Eyring model:

𝜀̇𝑝 𝑙𝑎𝑠𝑡𝑖𝑐_𝐹𝑖𝑏𝑟𝑖𝑙 = 𝑘 𝐵 𝑇 ℎ 𝑒 -∆𝐻 𝐹𝑖𝑏𝑟𝑖𝑙 𝑘 𝐵 𝑇 ⁄ sinh ( 𝛼 𝐹𝑖𝑏𝑟𝑖𝑙 𝜎 𝑘 𝐵 𝑇 ) (11)
The activation barrier of this secondary plasticity should be higher than that of the initial plasticity since it always occurs at higher stress level.

Finally, the behavior of the crystalline network branch is elastic then plastic with a strong softening due to the progressive collapse of the crystalline network. However, the overall behavior of the material does not exhibit softening but hardening.

As a consequence, it means that the second branch (consisting of the macromolecular network)

and partially the hardening of the confined amorphous phase, takes over on the crystalline network providing a strong hardening that balances the loss of stress due to the crystalline branch. Note that a similar idea was qualitatively proposed by Peterlin [START_REF] Peterlin | Plastic deformation of polymers with fibrous structure[END_REF].

Numerical identification of the model parameters

Following the physical considerations described above, we finally need to solve a two-body system explicitly defined by the following equations:

𝜎 = 𝜎 𝑐𝑟𝑦𝑠𝑡𝑎𝑙 𝑠𝑡𝑎𝑐𝑘 𝑛𝑒𝑡𝑤𝑜𝑟𝑘 + 𝜎 𝑚𝑎𝑐𝑟𝑜𝑚𝑜𝑙𝑒𝑐𝑢𝑙𝑎𝑟 𝑛𝑒𝑡𝑤𝑜𝑟𝑘 (12)
with 𝜎 𝑐𝑟𝑦𝑠𝑡𝑎𝑙 𝑠𝑡𝑎𝑐𝑘 𝑛𝑒𝑡𝑤𝑜𝑟𝑘 defined by equations 7-11 and 𝜎 𝑚𝑎𝑐𝑟𝑜𝑚𝑜𝑙𝑒𝑐𝑢𝑙𝑎𝑟 𝑛𝑒𝑡𝑤𝑜𝑟𝑘 defined by equations 2-6.

Numerically, the time integration of variables ePl (eq.9), ePl_Fibril (eq.11) and eflow (eq.6) is performed using an order two Runge-Kutta scheme.

The overall deformation of the semicrystalline polymer could be separated in several steps associated with different physical parameters such as crystallinity and entanglement network. The numerical identification of the parameters has been thus separated in three steps.

To take into account the dispersion of experimental results, the parameter identification was made with target values corresponding to the stress at selected deformations: eight values per strain rate.

Each of these values is an average obtained from the four compressive strain-stress curves for each strain rate.

In the first step, the Young's modulus (describing the crystal network contribution) has been determined by a large set of polyethylene with various crystallinity. Here we use the combination of experimental values from Humbert et al. [START_REF] Humbert | A re-examination of the elastic modulus dependence on crystallinity in semi-crystalline polymers[END_REF] and our own values (averaged over four experiments) which give access to parameters 𝑘 𝑖𝑛𝑖𝑡 (eq.10), 𝜇 𝐶𝐴 and 𝑁 𝐶𝐴 (eq. 8).

Parameters 𝜇 𝐶𝐴 and 𝑁 𝐶𝐴 have been verified on the unload part of the compressive curves (immediate recovery).

In the second step, the parameters of the plastic deformation (∆𝐻 𝑃𝑙 𝑎𝑛𝑑 𝛼 𝑃𝑙 , in equation 9) have been identified on the evolution of the yield point as a function of the strain rate for the weight 10.5 Mg/mol. Here also, the average values based on four experimental tests are used.

In the third step, the strain-hardening of the semicrystalline polymer is determined by a combination of three phenomena:

-The decrease of the intra-cohesive lamellar forces of the crystal network due to the shear of the crystal lamellae and the destruction of the spherulitic structure in favor of a fibrillar structure. In the model proposed here the parameter 𝜀 𝑁𝑎𝑡_𝐶𝑟𝑦𝑠𝑡𝑎𝑙 (eq.10), corresponding to the deformation for which the fibrillary state is reached, allows to reproduce this decrease.

-The strain hardening due to the entanglement network modeled here by the entropic elasticity decomposed in two parts: the elastic modulus 𝐶 0 (parameters 𝑁 𝑀𝑁 and 𝜇 𝑀𝑁 of the Langevin spring in equation 2) and the dependence with strain rates (parameters ∆𝐻 𝑀𝑁 𝑎𝑛𝑑 𝛼 𝑀𝑁 that appears in the expression of 𝐶 𝑎𝑝𝑝 (𝜀)in equation 3) combine with the effect of the confined amorphous phase.

-And finally, the disentanglement of the macromolecular chains described by the reptation theory (𝜂 𝑓𝑙𝑜𝑤 = 𝐺 𝑓𝑙𝑜𝑤 𝑀 3 ) and the parameters ∆𝐻 𝐹𝑖𝑏𝑟𝑖𝑙𝑠 and 𝛼 𝐹𝑖𝑏𝑟𝑖𝑙𝑠 (eq. 11) describing the slipping of the fibrils.

Consequently, the 7 parameters ( 𝜀 𝑁𝑎𝑡_𝐶𝑟 , 𝑁 𝑀𝑁 , 𝜇 𝑀𝑁 , ∆𝐻 𝑀𝑁 𝑎𝑛𝑑 𝛼 𝑀𝑁 , ∆𝐻 𝐹𝑖𝑏𝑟𝑖𝑙𝑠 and 𝛼 𝐹𝑖𝑏𝑖𝑟𝑙𝑠 𝑎𝑛𝑑 𝐺 𝑓𝑙𝑜𝑤 ) have been identified on the evolution of the strain hardening with strain rates and for the molecular weight 10.5 Million g/mol. To obtain the fitting curves of the compressive tests on the material with molecular weights 0.6 and 3.9 Million g/mol only parameters 𝜀 𝑁𝑎𝑡_𝐶𝑟 , ∆𝐻 𝑀𝑁 and 𝑁 𝑀𝑁 have been numerically identified. All the other parameters are the same as in the case of the molecular mass 10.5 Million g/mol. The complete set of parameters obtained after optimization is given in Table 1.

Several numerical checks were made to verify the consistency of the model: distribution of the strain and stress in the various constituent parts and energy conservation in the system which was verified to be valid in the limit of numerical errors.

Discussion of the influence of molecular weight and crystallinity on the model parameters

Figure 8 shows the confrontation of the optimized numerical modeling and experimental data of compressive tests for the 3 molecular weights 0.6; 3.9; and 10.5 Million g/mol. This model permits to capture the increase of the 'strain-hardening modulus' with the strain-rate and the evolution of the recovery for the three molecular weights. Several verifications of the behavior of the model were made to confirm its consistency with the literature on the mechanical behavior of semicrystalline polymers. For instance, Figure 9 shows the evolution of the yield stress and the apparent natural deformation (corresponding to the deformation for which the fibrillar structure is reached) with the crystallinity. As expected, these two physical parameters increase with the crystallinity. This can be explained by the percolation of the crystal network which improves as the crystallinity increases. Figure 11 shows the results of an extrapolation of the model to a deformation of 80% for the molecular weights 0.6, 3.9 and 10.5 Mg/mol and strain rates between 10 -4 s -1 to 1000s -1 .

Here again, the results are in agreement with the experimental observations similarly to what was said regarding Figure 10 for the strain rate 0.1s -1 . Figure 13 shows an extrapolation of the model to larger deformation and to cyclic loading conditions. These tensile and compressive curves reveal a larger permanent deformation for the lower molecular weight. In other words, the lowest molecular weight accumulates plasticity with the cycles that could potentially lead to failure whereas the cyclic behavior of the highest molecular weight specimen is mostly reversible. This result is in agreement with the differences in cavitation erosion resistance observed for different molecular weights of the UHMWPE (Deplancke et al., 2015a). This promising typical feature allows to envisage this model for numerical modeling of cavitation erosion where repetitive impacts of multiple bubble collapses need to be taken into account.

CONCLUSION

An 1D hyperelastic-viscoplastic constitutive model is proposed to describe the mechanical compressive behavior of semicrystalline polymer with ultra-high molecular weight under a wide range of strain rates. The mechanical behavior of sintered UHMWPE with different molecular weight (0.6; 3.9 and 10.5 Mg/mol) was investigated over a large range of strain rates up to a deformation of 0.4. Experimentally, this range of strain was achieved using uniaxial compression tests performed on a servo-hydraulic testing machine for the lowest strain rates (10 -4 to 10 s -1 )

The mechanical behavior is described by a physical approach of the semicrystalline polymers based on a macromolecular network bridled by the crystalline network formed by lamellar stacking. In the case of UHMWPE, the impact of the macromolecular network, both long entangled chains as well as tie molecules, leads to a plastic deformation of the specimen without necking.

Indeed, the modulus collapse due to the splintering of the crystalline networks is immediately caught up by the active chains of the macromolecular network.

The originality of the proposed model lies in i) allowing the description of the crystalline network collapse by the decrease along the plastic deformation of the mechanical coupling between rubbery and crystalline phase i.e. the stress transfer to the hard/crystalline phase decreases as the plastic deformation increases, and ii) describing the role of both the confined amorphous phase and the long-entangled chains allowing strain hardening, short and long-term deformation recovery.

After a numerical identification of the twelve parameters of the model, a good agreement was observed between model response and experimental results in a large range of strain rates.

Furthermore, the same set of parameters allows to describe 3 different molecular weights (0.6, 3.9 and 10.5 Million g/mol), except two parameters which are logically different since the crystallinity of these materials is different: the percolation of the crystal network and the crystallinity. In addition, due to its physical consistency, the model gives satisfactory qualitative responses both for tensile and cyclic loading conditions largely above the strain range for which the parameters have been adjusted. This confirms that the essential of the physic required to describe the outstanding resistance to cavitation of the UHMWPE is captured by this simple 1D model.

Ongoing work consists of two main actions: (i) investigating other loading conditions such as cyclic loads and (ii) implementation of a 3D version of the model in the finite element code ABAQUS.

Figure 2 :

 2 Figure 2: Strain-stress curves of the experimental compressive tests for molecular weights 0.6 and 10.5 Mg/mol at strain rate 0.1s -1 .

Figure 3 .

 3 Figure 3. Schematic view of (a) the entanglement of crystals and macromolecular network and (b)

Figure 4 .

 4 Figure 4. Schematic view of the evolution of rheological model for large strain (fibrillation and stretched

(

  see fig 5). A viscous dashpot is added in series to the macromolecular network description in order to allow for possible creep of the macromolecular network since the chains are not cross-linked, and ultimately to describe the damage by viscous flow. The mechanism associated with the flow could be chains fractures or reptation activated by the stress.

Figure 5 .

 5 Figure 5. Rheological representation of the part "Macromolecular Network" of the constitutive

Figure 6 .

 6 Figure 6. Schematic view of the rheological representation of the part "Crystalline stacks

Figure 7 .

 7 Figure 7. Rheological representation of the part "Crystalline stacks network" of the constitutive

Figure 8 :

 8 Figure 8: Experimental and numerical strain-stress curvesfor compression test performed on the three molecular weights for various strain rates in the range from 10 -4 s -1 to 1000s -1 . The

  first two columns in bold italics are the molecular weights given by the manufacturer, Celanese, and the crystallinity measured by DSC. The two columns in red are the only parameters that change with the molecular weight.

Figure 9 :

 9 Figure 9: Evolution of the yield stress and the apparent natural deformation as a function of the

Figure 10 :

 10 Figure 10: Comparaison of experimental and numerical stress-strain curves obtained for compression tests performed at 0.1 s -1 on UHMPE with a molecular weight of (a) 0.6Mg/mol and

Figure 11 :

 11 Figure 11: Stress-strain curves of the numerical extrapolation of the model to a deformation of 80% for the three molecular weights and various strain rates in the range from 10 -4 s -1 to 1000s -1 .

Figure 12 :

 12 Figure 12: Tensile experiments (black dots) and extrapolation of the model (red curves) to

Figure 13 .

 13 Figure 13. Extrapolation of the model for cyclic loading (10 cycles) (a) in compression for different strain rates and (b) in tensile for a strain rate 0.001 s -1 .

Table 1 :

 1 Optimized set of parameters of the constitutive model for the crystal network and for the confined amorphous phase (first table) and to the macromolecular network (second table). The

complemented by magnesium split Hopkinson pressure bar tests for the highest strain rates (10 3 s - 1 ).

ACKNOWLEDGEMENTS

This research was carried out under project CARAPACE ANR-16-ASTR-0004-04 which provided financial support for Dr. T. Deplancke.

The authors are indebted to Celanese (Oberhausen, Germany) for the generous supply of the UHMWPE sample together with its molecular characteristics. The authors also gratefully acknowledge the support of N. Ranc for the realization of the SHPB tests.


Carmen Stelian
email: steliancarmen@yahoo.com

Jessica Muzy

Serge Labor

Marc Fivel

Hugues Cabane

Thierry Duffar

Numerical Analysis of Thermal Stress in Semi-Transparent Oxide Crystals Grown by Czochralski and EFG Methods

Keywords: Czochralski, single crystal growth, stresses, sapphire, computer simulation

Introduction

Oxide crystals such as sapphire or langatate (La3Ta0.5Ga5.5O14) are used as raw materials for acustooptic-electronic devices. These oxides are usually grown by the Czochralski (Cz), Edge-defined Film-fed Growth (EFG), Heat Exchanger Method (HEM) or Kyropoulos techniques. Crystal cracking and high dislocation density are often reported in Cz growth of oxides [1][2][3][4][5]. Numerical investigations of thermal stresses in Cz growth of oxide crystals have shown that crystal cracking and dislocations are related to the temperature field in the furnace and the shape of the crystal-melt interface [6][7][8]. Oxide crystals are semi-transparent to the infrared radiation by absorbing and emitting the radiative heat. Therefore, the temperature field in the crystal and the shape of the growth interface are significantly affected by the internal radiative heat transfer. Deep convex shapes of the crystal-melt interface are experimentally and numerically reported in these systems [9][10][11][12]. Crystal rotation is usually applied to flatten the growth interface and to reduce thermal stresses [11][12][13][14]. Recently, Stelian et al. [11,12] revealed by means of numerical modeling that applying crystal/crucible rotation in some systems results in distorted shapes of the growth interface and increased stresses in the ingots.

The effect of anisotropic and temperature dependent mechanical properties on thermal stresses distribution in oxide crystals was investigated in [10,12]. It was found that the anisotropic effect is significant in LiNbO3 crystals [10] but less important in the sapphire crystals [12].

The objective of the present work is to compare the effect of anisotropic mechanical properties on thermal Crystal growth of oxides is generally difficult since large curvatures of the growth interface in these systems generate high thermal stress, dislocations and crystal cracking. Three-dimensional numerical modeling is applied to investigate thermal stress distribution in sapphire and langatate La3Ta0.5Ga5.5O14 (LGT) semi-transparent crystals grown by Czochralski (Cz) and Edge-defined Film-fed Growth (EFG) techniques. The analysis of thermal stress distribution in a sapphire ingot of 5 cm in diameter grown in a Czochralski furnace shows high von Mises stress distributed on large areas in the internal part of the crystal. The distribution of thermal stress is almost axisymmetric, which means that the effect of the anisotropic mechanical properties of sapphire is low. Thermal stress computations for piezoelectric langatate crystals of 5 cm in diameter grown in a Czochralski configuration show non-symmetrical von Mises distribution with higher stress on one side of the ingot. These numerical results are in agreement with experimental results showing non-symmetrical cracking at the outer surface of the crystal. 3D modeling is applied to investigate thermal stress in a multi-die EFG growth system of white sapphire ribbons. Computations show that the von Mises stress is almost constant with increasing the number of ribbons from 2 to 10. Two models are applied to simulate the internal radiative heat transfer in the sapphire crystals: P1 approximation and the Rosseland radiation model. Numerical results show that applying Rosseland formula introduces significant errors in temperature field calculations especially in the case of the EFG configuration.

stresses distribution in two oxide crystals (sapphire and LGT) grown by the Czochralski method. Thermal stresses will be also investigated in the case of sapphire ribbons grown in a multi-die EFG system. Finally, the application of two models (P1 approximation and Rosseland radiation model) to simulate the internal radiative heat transfer in sapphire crystals will be discussed.

Model description

The coupled equations describing the heat transfer, momentum and thermal stresses are solved by using the finite element software COMSOL Multiphysics:

𝜌𝑐 𝑃 ( 𝜕𝑇 𝜕𝑡 + 𝑢 ⃗ • ∇𝑇) = 𝑘∇ 2 𝑇 + 𝑄 𝑅 (1) 𝜌[ 𝜕𝑢 ⃗ ⃗ 𝜕𝑡 + (𝑢 ⃗ • ∇)𝑢 ⃗ ] = -∇𝑝 + 𝜇∇ 2 𝑢 ⃗ -𝜌𝑔 [𝛽 𝑇 (𝑇 -𝑇 𝑟𝑒𝑓 )](2) 𝜎 𝑖𝑗 = 𝐶 𝑖𝑗𝑘𝑙 [𝜀 𝑘𝑙 -𝛼 𝑘𝑙 (𝑇 -𝑇 𝑟𝑒𝑓 ] (3)
In the above equations, 𝜌 , 𝑐 𝑝 , 𝑘 , 𝑄 𝑅 , 𝑝 , 𝜇 , 𝛽 𝑇 , 𝑇 , 𝑢 ⃗ , 𝑔 , 𝜎 𝑖𝑗 , 𝜀 𝑘𝑙 , 𝐶 𝑖𝑗𝑘𝑙 , 𝛼 𝑘𝑙 and 𝑇 𝑟𝑒𝑓 denotes the density, specific heat, thermal conductivity, radiative heat source, pressure, dynamic viscosity, thermal expansion coefficient, temperature, flow velocity, gravity vector, stress tensor components, strain tensor components, components of the elasticity tensor, components of thermal expansion tensor and the reference temperature.

In the case of Czochralski method, temperatures carried out from 2D-axisymmetric global modeling of the furnace are used as boundary conditions for three-dimensional (3D) local modeling of the stress distribution in the crystal. The EFG multi-die system has two orthogonal symmetry planes, so, in this case, 3D global modeling is applied to compute the temperature field in one fourth of the furnace. Then, 3D local modeling is applied to compute the stress distribution in the sapphire ribbons.

A structured/unstructured combined scheme is applied in 2D global modeling of the furnace. The triangular mesh consists of maximum 15000 domain elements. The mesh used in 3D modeling is unstructured and contains up to 150000 tetrahedral elements. The minimum element size used in 3D thermal stress computations is 0.005 cm.

Temperature computations include heat transfer by conduction, wall-to-wall radiative exchanges and internal radiative heat transfer in the semi-transparent crystals. Two models are applied (and compared) to simulate the internal radiation effect in sapphire crystals grown by Cz and EFG methods: P1 approximation and Rosseland radiation model.

P1 approximation model considers a radiative heat source 𝑄 𝑅 in Eq. ( 1), which is given by:

𝑄 𝑅 = 𝑎(𝐺 -4𝑛 2 𝜎𝑇 4 ) (4)
In the above equation, 𝑎 is the absorption coefficient, 𝐺 is the incident radiation, 𝑛 is the refractive index and 𝜎 is the Stefan-Boltzmann constant. The melt was considered opaque for both sapphire and LGT computations. The Rosseland radiation model, which is frequently used in numerical modeling of oxides crystals growth [15,16], assumes an increased thermal conductivity of the semi-transparent crystal:

𝑘 𝑒𝑓𝑓 = 𝑘 𝑜𝑝𝑎𝑞𝑢𝑒 + 16𝑛 2 𝜎𝑇 3 3𝑎 (5)
𝑘 𝑒𝑓𝑓 is the effective thermal conductivity of the semitransparent crystal and 𝑘 𝑜𝑝𝑎𝑞𝑢𝑒 is the thermal conductivity of the opaque crystal.

Heat transfer computations include the latent heat released at the growth interface and the effect of the crystal translation. The shape of the crystal-melt interface is computed by using a deformable mesh technique only in 2D-axisymmetric computations. In the case of 3D modeling of the EFG multi-die system, the shape of the crystal-melt interface was considered flat.

Flow computations include buoyancy convection, Marangoni effect and the forced convection generated by the crystal rotation. The flow is laminar for both Cz simulated cases. The Grashof numbers (𝐺𝑟 = 𝑔𝛽 𝑇 ∇𝑇𝑅 3 𝜌 2 /𝜇) are: 𝐺𝑟 = 1200 (for sapphire growth) and 𝐺𝑟 = 4000 (for LGT growth).

Thermal stress calculations are three-dimensional and consider the effect of the anisotropic elastic constants. Sapphire is classified as a trigonal class 3𝑚 ̅ crystal with six independent elastic constants: 𝐶 11 = 497 𝐺𝑃𝑎, 𝐶 12 = 162.66 𝐺𝑃𝑎, 𝐶 13 = 117.26 𝐺𝑃𝑎, 𝐶 14 = 22.96 𝐺𝑃𝑎 , 𝐶 33 = 501.8 𝐺𝑃𝑎 , 𝐶 44 = 147.2 𝐺𝑃𝑎 [2].

LGT is a class 32 crystal and has six independent elastic constants: 𝐶 11 = 192 𝐺𝑃𝑎 , 𝐶 12 = 111 𝐺𝑃𝑎 , 𝐶 13 = 102.9 𝐺𝑃𝑎 , 𝐶 14 = 13.8 𝐺𝑃𝑎 , 𝐶 33 = 264.1 𝐺𝑃𝑎 , 𝐶 44 = 51 𝐺𝑃𝑎 [START_REF] Davulis | Characterization of the Elastic, Piezoelectric, and Dielectric Properties of Langatate At High Temperatures Up To 900 •C[END_REF]. In the present simulations, the z-axis is parallel to the threefold c-axis, and x-axis is parallel to the twofold a-axis of the crystals. The effect of temperature dependent elastic constants and the anisotropy in the thermal expansion coefficient was discussed in a previous work [12]. More details about the numerical and physical model used to simulate crystal growth of oxides are given in [11,12]. The physical properties Copyright line will be provided by the publisher of sapphire and LGT are taken from [2,[START_REF] Chen | [END_REF], and respectively [START_REF] Davulis | Characterization of the Elastic, Piezoelectric, and Dielectric Properties of Langatate At High Temperatures Up To 900 •C[END_REF]19].

Results and discussion

Crystals grown by the Czochralski technique

Stress distribution in the sapphire crystals

The modeling is performed for sapphire crystals of 5 𝑐𝑚 in diameter grown along c-axis in an inductive Czochralski furnace. Fig. 1a shows the temperature and the velocity fields carried out from 2D-axisymmetric simulations performed at a given stage of the growth process (the length of the cylindrical part of the crystal is 𝐿 𝑐𝑦𝑙 = 1.6 𝑐𝑚). The shape of the crystal-melt interface was computed by using the deformable mesh technique. Numerical results reveal a deep convex shape of the interface. The computed dimensionless interface deflection 𝜏 = 𝑓 𝑅 𝐶 = 3.6 𝑐𝑚 2.5 𝑐𝑚 = 1.44 ⁄ ⁄ (with 𝑓interface deflection and 𝑅 𝐶 -crystal radius) is in agreement with the experimental result 𝜏 𝑒𝑥𝑝 = 𝑓 𝑅 𝐶 = 6 𝑐𝑚 3.8 𝑐𝑚 = 1.57 ⁄ ⁄ of Bruni [3]. The curved shape of the isotherms in the crystal indicates large radial temperature variations that generate high thermal stresses in the sapphire ingot. The flow pattern exhibits one single vortex, with a maximum velocity of 𝑢 𝑚𝑎𝑥 = 0.03 𝑚/𝑠 at the symmetry axis. The isotherms are significantly distorted by the melt flow.

Fig. 1b shows the von Mises stress distribution at the outer surface of the crystal. The results are carried out from 3D thermal stress computations that use the temperatures carried out from 2D-axisymmetric global modeling of the furnace. High thermal stresses are observed at the triple crystal-melt-gas line and at the outlet surface of the conical part of the crystal.

Fig. 2 shows the contour lines of von Mises stress distribution plotted in the longitudinal section (x-z plane) and the transversal section (x-y plane at 𝑧 = 7 𝑐𝑚) of the crystal. The contour lines of the von Mises stress are almost axisymmetric in the longitudinal section (Fig. 2a). Large thermal stresses (greater than 200 𝑀𝑃𝑎) are observed in almost whole crystal volume. The maximum stress value is located at the triple crystal-melt-gas line (900 𝑀𝑃𝑎), while the minimum value (2.3 𝑀𝑃𝑎) is located in the seed. The stress distribution in the transversal section of the crystal reveals a slight asymmetry (Fig. 2b). The difference between the maximum and the minimum values of the von Mises stress at the outer surface of the crystal, Δ𝜎 𝑀𝑖𝑠𝑒𝑠 = 𝜎 𝑀𝑖𝑠𝑒𝑠 𝑚𝑎𝑥 -𝜎 𝑀𝑖𝑠𝑒𝑠 𝑚𝑖𝑛 = 31 𝑀𝑃𝑎 represents only 7% from the maximum stress value. The radial stress distribution along the x-axis in the transversal section of the crystal at 𝑧 = 7 𝑐𝑚 is plotted in Fig. 2c. High von Mises stress (over 160 𝑀𝑃𝑎) are distributed in a region of 5 𝑚𝑚 thickness located near the outer surface of the crystal. The critical resolved shear stress varies between 5 𝑀𝑃𝑎 and 110 𝑀𝑃𝑎 over the temperature range computed in the crystal [20]. That explains high dislocation densities measured in the sapphire crystals grown by the Czochralski technique. Fig. 5a shows that the maximum stress value is located at the outer surface of the crystal at ≈ 0.25 𝑐𝑚 distance above the triple crystal-melt-gas line in the direction parallel to y-axis. The minimum stress value at the external surface of the crystal is observed along the x-axis. Fig. 5b shows a LGT ingot of 5.4 𝑐𝑚 in diameter, which was decanted from the melt in order to visualize the shape of the growth interface. The crystal exhibits non-symmetrical surface cracking.

Sapphire crystals grown by the EFG technique

Sapphire ribbons of dimensions 0.8 𝑐𝑚 × 4 𝑐𝑚 × 14.5 𝑐𝑚 are grown in a multi-die EFG system. 3D global modeling is firstly applied to compute the temperature field in the furnace. To minimalize the computational time, the modeling is performed only for a quarter of the real domain by neglecting the convection in the melt. The shape of the crystal-melt interface is not computed in this case, since 3D computations with deformable mesh are too complicated. The meniscus height ℎ 𝑚 = 0.07 𝑐𝑚 was estimated by applying the analytical formula proposed in [21]. Fig. 6 shows the temperature field carried out from the global modeling of an EFG system with ten molybdenum dies. It can be observed a significant amount of heat lost by the crystals due to the internal radiative effect. The temperatures decrease by Δ𝑇 = 550 𝐾 from the bottom part to the top part of the crystals. Even if the vertical temperature gradients are high in this configuration, the horizontal temperature gradients (which mainly generate the thermal stresses) are minimized due to the heat radiative exchanges between the sapphire ribbons, that homogenize the temperatures in the horizontal direction.

Temperatures computed for EFG systems with 2, 4, 6, 8 and 10 dies are used in the local simulation of the thermal stresses distribution in the sapphire ribbons. The results plotted in Fig. 7 show that the maximum von Mises stress (≈ 60 𝑀𝑃𝑎) remains almost constant with increasing the number of ribbons from 2 to 10. The von Mises stress is higher in the vicinity of the crystal-melt interface, and slightly increases from the center to the outer surface of the crystal. The present results are comparable to the numerical computations of Théodore et al. [22] showing the same von Mises stress distribution map and a maximum value of ≈ 16 𝑀𝑃𝑎 in the case of a smaller sapphire ribbon.

Maximum stress values computed in the EFG multi-die system (≈ 60 𝑀𝑃𝑎) are smaller compared to results obtained in the case of Cz sapphire ingot, which revealed high thermal stresses (over 200 𝑀𝑃𝑎) spread in almost whole volume of the crystal.

Application of P1 approximation and Rosseland model to simulate the internal radiative heat transfer

Two models are usually applied to simulate the internal radiative heat transfer in semi-transparent crystals: P1 approximation and Rosseland radiation model. The Rosseland model, which is most frequently used to simulate the internal radiative effect in sapphire crystals, is based on the artificially enhanced thermal conductivity of the crystal, given by Eq. 5. This model was used for numerical modeling of heat transfer in participating media in [15,16] but the results for the curvature of the growth interface show significant deviations from the experimental observations. The P1 approximation model takes into account an additional heat source term in Eq. 1 that describes the radiative heat transfer for an absorbing, emitting and scattering medium (the scattering effect is not considered in the present simulations). Numerical modeling based on this model has shown that deep convex shapes of the crystal-melt interface are explained by the internal radiative heat exchanges [9,11,23,24]. Numerical simulations in [11,23,24] have been validated by comparing the computed shape of the growth interface to the experimental results.

In the present section, P1 approximation and the Rosseland model are applied to simulate the heat transfer in Cz and EFG growth of sapphire crystals. Fig. 8 compares the numerical results carried out from Cz simulations which use the Rosseland model (Fig. 8a) and the P1 approximation (Fig. 8b). The effective thermal conductivity of the sapphire crystal (𝑘 𝑒𝑓𝑓 ≈ 500 𝑊/𝐾 • 𝑚) was estimated by using Eq. 5 (with 𝑎 = 20 𝑚 -1 and 𝑛 = 1.77). The shape of the growth interface was computed by using the deformable mesh technique. During a transient time (50 -300 𝑠), the mesh is deformed in such way as the surface between the crystal and the melt (shown on the figure at time = 0 𝑠 ) takes the shape of the solidification isotherm. Fig. 8a shows the numerical plot at the beginning ( 𝑡 = 0 𝑠) and the end ( 𝑡 = 300 𝑠) of the transient computation with deformable mesh. The dimensionless interface deflection 𝜏 = 𝑓 𝑅 𝐶 = 2.4 𝑐𝑚 2.5 𝑐𝑚 = 0.96 ⁄ ⁄ is smaller than observed in experiments [3]. Temperatures in the crystal are homogeneous and close to the melting point, which indicate that the amount of heat lost by the crystal is very small.

The results carried out by applying the P1 approximation show much curved interface (Fig. 8b). The dimensionless interface deflection 𝜏 = 𝑓 𝑅 𝐶 = 3.6 𝑐𝑚 2.5 𝑐𝑚 = 1.44 ⁄ ⁄

, is much closer to the experimental result of Bruni [3]. The temperatures in the crystal decrease significantly due to the internal radiative effect. That generates large temperature gradients in the crystal, and, consequently, high thermal stresses, dislocations and cracking.

Numerical results obtained by applying the Rosseland and P1 approximation models in the case of EFG growth of 10 sapphire ribbons are plotted in Fig. 9a and, respectively, Fig. 9b. The crystal-melt interface was considered plane in these simulations. The power in the furnace was varied until the melting point is reached in the liquid meniscus corresponding to the die 6, which is nearest from the symmetry plane (central dies are labeled as 5 and 6). Fig. 9a shows the temperature distribution in the vertical (x-z) section at the center of the dies (symmetry plane). The red line represents the solidification isotherm. It can be observed that the temperatures decrease significantly in the horizontal direction (from the ribbon 10 to the ribbon 6). The solidification isotherm is located in the liquid meniscus only for the dies 6-8, then deviates far from the liquid meniscus for the dies 9 and 10. That should make impossible the simultaneous growth of 10 ribbons, which is in contradiction with our experimental results. Fig. 9b shows more realistic results carried out from simulations which use the P1 approximation. In this case, the temperatures decrease on the vertical direction, from the bottom part to the top part of the ribbons. The solidification isotherm is located near the crystal-melt interfaces for all the ribbons.

Conclusions

Three-dimensional numerical simulations were performed to investigate the effect of the anisotropic elastic constant on thermal stresses distribution in sapphire and piezoelectric langatate crystals grown by the Chochralski and EFG techniques.

The modeling performed in the case of sapphire crystals grown in an inductive Cz configuration shows large thermal stresses distributed almost symmetrically in the whole volume of the crystal. This result could explain high dislocation densities measured in the sapphire crystals grown by this technique.

Copyright line will be provided by the publisher

The computations performed in the case of langatate crystals grown in a Cz configuration show non-symmetrical stress distribution with maximum stress located at one side of the crystal. The von Mises stresses are relatively small in the internal part of the crystal, but increase significantly in a thin region of 1-2 mm located near the outer surface of the crystal. This explains non-symmetrical cracks experimentally observed at the external surface of the crystal. Large temperature gradients at the outer surface of the crystal could be generated by the additional thermal shields inserted around the crystal in this configuration.

Thermal stress computations performed in the case of a multi-die EFG system show that the distribution and the maximum von Mises values are almost unchanged with increasing number of ribbons from 2 to 10. This is explained by the radiative heat exchanges between the ribbons that homogenize the temperatures in the horizontal direction. The maximum stress values are smaller as compared to the case of Cz sapphire growth, which predicts lower dislocation densities in sapphire ribbons grown by the EFG technique.

The application of two models (P1 approximation and Rosseland model) to simulate the internal radiative heat transfer was investigated in the case of sapphire crystals grown by Cz and EFG techniques. The present results show that the Rosseland model introduces significant errors in heat transfer computations, especially in the case of the EFG system. It is found that the interface curvatures are generally underestimated when the Rosseland model is applied to simulate Cz growth of sapphire. EFG thermal computations using the Rosseland model give unrealistic results for temperatures distribution and the position of the solidification isotherm. On the other hand, the application of P1 approximation model give more realistic results, which are validated by comparison to experiments. It is concluded that the Rosseland model, which is still frequently used in numerical works, is not appropriate to simulate the internal radiative effect in semi-transparent crystals. Original Papers

Table of Contents

Figure 1

 1 Figure 1 Numerical results obtained in the case of Czochralski sapphire growth: (a) temperature distribution and isotherms (left side), velocity field (right side) carried out from global 2D axisymmetric computations. (b) von Mises stress distribution at the surface of the crystal carried out from local 3D computations.

Figure 2

 2 Figure 2 Numerical results obtained in the case of Czochralski sapphire growth: (a) von Mises stress distribution in the longitudinal (xz) section of the crystal (maximum and minimum values are shown on the figure); (b) von Mises stress distribution in the transversal (xy) section of the crystal at 𝑧 = 7 𝑐𝑚 (maximum and minimum values are shown at the outer surface of the crystal); (c) von Mises stress distribution along the x-axis in the transversal (x-y) section of the crystal at 𝑧 = 7 𝑐𝑚.

Figure 3

 3 Figure 3 Numerical results obtained in the case of Czochralski LGT growth: (a) temperature distribution and isotherms (left side), velocity field (right side) carried out from global 2D axisymmetric computations. (b) von Mises stress distribution at the surface of the crystal carried out from local 3D computations.

Figure 4

 4 Figure 4 Numerical results obtained in the case of Czochralski LGT growth: (a) von Mises stress distribution in the longitudinal (x-z) section at 𝑦 = 1𝑐𝑚 (maximum and minimum values are shown on the figure); (b) von Mises stress distribution in the transversal (x-y) section of the crystal at 𝑧 = 6.6 𝑐𝑚 (maximum and minimum values are shown at the outer surface of the crystal); (c) von Mises stress distribution along x-axis in the transversal (x-y) section of the crystal at 𝑧 = 6.6 𝑐𝑚.

Figure 5

 5 Figure 5 (a) Von Mises stress distribution at the outer surface of the LGT crystal. Maximum and minimum stress values are plotted in a transversal section at 0.25 cm distance above the crystalmelt-gas line; (b) LGT ingot decanted from the melt which exhibits non-symmetrical surface cracking.

Figure 6

 6 Figure 6 Temperature field carried out from the global modeling of an EFG system with 10 molybdenum dies. Only a quarter of the real domain was simulated.

Figure 7

 7 Figure 7 Von Mises stress distribution in sapphire ribbons grown in an EFG multi-die system with 2, 4, 6, 8 and 10 dies. Red color indicates 50 MPa.

Figure 9

 9 Figure 9 Temperatures computed in the case of EFG growth of 10 sapphire ribbons by applying: (a) Rosseland radiation model; (b) P1 approximation. The red line represents the solidification isotherm.

  be provided by the publisher

  Page xx-yy (for typesetter) C. Stelian, J. Muzy, S. Labor, M. Fivel, H. Cabane, and T. Duffar Numerical analysis of thermal stress in semi-transparent oxide crystals grown by Czochralski and EFG methods Modeling of Czochralski growth of sapphire crystals shows large thermal stresses distributed symmetrically in the whole volume of the crystal. That could explain high dislocation densities measured in the sapphire ingots. Computations performed in the case of langatate crystals show non-symmetrical stress distribution, with maximum stress located at one side of the crystal. That explains non-symmetrical surface cracking observed in the experiments. Numerical results show that the Rosseland model is not appropriate to simulate the internal radiative effect in semi-transparent crystals.

Copyright line will be provided by the publisher