J Andrieux

O Dezellus

F Bosselet

J C Viala

Low-Temperature Interface Reaction between Titanium and the Eutectic Silver-Copper Brazing Alloy

Keywords: metals, intermetallics, diffusion paths, chemical potential gradients, ternary system

Reaction zones formed at 790 °C between solid titanium and liquid Ag-Cu eutectic alloys (pure and Ti-saturated) have been characterized. When pure Ag-Cu eutectic alloy with 40 at.% Cu is used, the interface reaction layer sequence is: αTi / Ti 2 Cu / TiCu / Ti 3 Cu 4 / TiCu 4 / L. Because of the fast dissolution rate of Ti in the alloy, the reaction zone remains very thin (3-6 µm) whatever the reaction time. When the Ag-Cu eutectic alloy is saturated in titanium, dissolution no longer proceeds and a thicker reaction zone with a more complex layer sequence grows as the reaction time increases. Four elementary chemical interaction processes have been identified in addition to Ti dissolution in the liquid alloy. These are growth of reaction layers on Ti by solid state diffusion, nucleation and growth from the liquid of TiCu 4 , isothermal solidification of silver and, finally, chemical conversion of the Cu-Ti compounds by reaction-diffusion in the solid state. A mechanism combining these processes is proposed to account for the constitution of Ti/ Ag-Cu/ Ti joints brazed at 780-800 °C.

Introduction

The Ag-Cu eutectic alloy with 40 at.% Cu (28 wt.% Cu ) has already shown good wetting and bonding characteristics with regard to conventional brazing of titanium parts at temperatures higher than 850 °C [START_REF] Eustathopoulos | Wettability at High Temperature[END_REF]. This behaviour is explained by the formation of binary Ti-Cu compounds that promote reactive wetting and spreading at the liquid/solid interface [START_REF] Eustathopoulos | Wettability at High Temperature[END_REF][START_REF] Wu | Growth Mechanism of Compounds at the Interface of Ti6Al4V/Ag-28Cu Filler Alloy[END_REF][START_REF] Shiue | Interfacial Reactions of IR Brazing Cu and Ti[END_REF][START_REF] Liu | The Microstructural Observation and Wettability Study of Brazing Ti-6Al-4V and 304 Stainless Steel Using Three braze Alloys[END_REF]. Surprisingly, data are very scarce about brazing of titanium parts at temperatures lower than 850 °C, although promising observations have been made [START_REF] Shiue | Interfacial Reactions of IR Brazing Cu and Ti[END_REF]. In fact, operating at temperatures closer to the melting point of the alloy (which is of 780 °C) should be favourable to the preservation of the microstructure and composition of the titanium base metal. Moreover, formation of copper-rich Ti-Cu solid compounds should render possible isothermal solidification of the Ag-Cu brazing alloy. In other words, the development of a process for brazing titanium or titanium alloys by transient liquid phase bonding at low temperature should be very attractive from an applied perspective [START_REF] Macdonald | Transient liquid Phase Bonding[END_REF][START_REF] Macdonald | Low Temperature Transient Liquid Phase Bonding of Ti-6Al-4V, International Trends in Welding Science and Technology[END_REF]. For these reasons, it was decided to further investigate the chemical interaction between titanium and the eutectic Ag-Cu alloy at temperatures down to 780 °C, the temperature of the eutectic transformation. In the first approach at investigating the kinetics of compound formation in the Ag-Cu-Ti system, cold-pressed mixtures of Ag, Cu and Ti powders were heated at 700-860 °C and characterized [START_REF] Andrieux | Details on the Formation of Ti 2 Cu 3 in the Ag-Cu-Ti System in the Temperature Range 790 to 860[END_REF]. In the present work, additional experiments have been carried out at 780-800 °C on Ti/Ag-Cu semi infinite couples and Ti/Ag-Cu/Ti sandwiches in order to determine the extent and constitution of the interface reaction zone formed between solid Ti and liquid Ag-Cu eutectic alloy.

Experimental

Experiments on Ti/Ag-Cu couples were performed in the sessile drop configuration. On the one hand, a pure binary Ag-Cu eutectic alloy (40 at.% Cu ) without Ti was used; on the other hand, 1.88 at.% of Ti (1 wt.%) was added so that the Ag-Cu alloy was slightly oversaturated in Ti at 790 °C according to [START_REF] Eremenko | The liquidus Surface of the System Ti-Cu-Ag[END_REF]. In both cases, the quantities used were such that the composition of the liquid remained as constant as possible. For that purpose, small blocks of Ag-Cu or Ag-Cu-Ti alloy (5 mm x 3 mm x 4-5.5 mm, 0.6-0.82 g) were prepared by coldpressing (270 MPa) commercial powders of silver (99,99 wt.% pure, Goodfellow), copper (99 wt.% pure, Goodfellow) and titanium (purity higher than 98 wt.%, grain size after sieving d < 100 µm, Merck). The blocks were put on a freshly polished and degreased titanium plate (10 mm x 5 mm x 2 mm) that was entirely coated with a Y 2 O 3 -based oxide liner (Stopyt, Morgan), except for a square area of 2 mm x 2 mm at its upper face that remained uncovered and was thus in direct contact with the Ag-Cu(-Ti) alloy. The samples were put in an alumina boat and placed inside a silica reaction tube evacuated under a dynamic primary vacuum of ∼1 Pa. Another boat filled with titanium powder was also placed in the hot zone of the reaction tube to act as a getter protecting the samples from oxidation, nitridation or carburization. The reaction tube was heated using a horizontal tubular furnace regulated to better than ± 1 °C. The temperature was first raised in less than 30 min from 25 °C to 790 °C and then held at that value for durations varying from 2 to 50 min. At the end of the isothermal heat-treatment, the samples were cooled to room temperature by pulling the silica reaction tube out of the furnace. The mean cooling rate during the first 100 °C drop was faster than 10 °C s -1 . The heat-treated samples were diamond sawn, polished and examined by optical metallography (OM) and field-emission gun scanning electron microscopy (FEG-SEM, Hitachi S-800). Analyses were carried out by electron probe microanalysis (EPMA, Cameca Camebax) with an accelerating voltage of 10 kV and a beam current of 9 nA: under these conditions the lateral resolution was of ±0.7 µm. After background subtraction, the counting rates obtained for Ag, Cu and Ti in at least ten different points of the same layer or type of crystals were averaged and referred to the counting rates recorded under the same conditions on pure and freshly polished element standards. After corrections for atomic number, absorption and fluorescence, the atomic contents of Ag, Cu and Ti were obtained with an accuracy better than ± 0.5 at%. It had been previously verified by X-ray powder diffraction [START_REF] Andrieux | Details on the Formation of Ti 2 Cu 3 in the Ag-Cu-Ti System in the Temperature Range 790 to 860[END_REF] that the crystalline phases for which each interface reaction layer or crystal analysed were effectively obtained under the same heat-treatment conditions.

Results

Depending on whether the Ag-Cu alloy was saturated in titanium or not, the reaction zones formed between solid titanium and the liquid alloy varied in thickness, morphology and composition.

1. Ag-40 at.% Cu liquid saturated in Ti

When titanium reacted at 790 °C with Ti-saturated Ag-Cu eutectic alloy, a thick reaction zone rapidly developed at the liquid/solid interface. After 2 min at 790 °C, the reaction zone had a total thickness of about 15 µm. After 15 and 50 min at that temperature, the reaction zone thickness attained ~100 µm and ~150 µm, respectively. For the latter two reaction times, the reaction zone morphology was the same. As shown in Fig. 1, five successive layers can be distinguished between the unconverted titanium substrate and the solidified Ag-Cu liquid saturated in Ti. In direct contact with titanium, a first thin layer (less than 2 µm thick) with a single-phased appearance has formed. According to its EPMA analysis (Table 1), this first layer consists of the phase Ti 2 Cu with some silver in solid solution. A second layer (less than 3 µm thick) surrounding the first one also exhibits a single-phased appearance: it analyzes by EPMA (Table 1) for the compound TiCu with some silver in solid solution. A third layer which is thicker than the two previous ones (about 10 µm thick) and which clearly consists of two different phases can then be observed: in that layer, the dark phase analyzes for Ti 3 Cu 4 whereas the bright inclusions are made of metallic silver with about 9 at% of copper and 5 at% of titanium dissolved (Table 1). The fourth layer rapidly contains more bright phase than the third so that the bright phase becomes a matrix which embeds dark elongated crystals. According to EPMA (Table 1), the bright phase consists of metallic silver with 9-14 at.% of copper and 5-1 at.% of titanium dissolved and the dark phase is made of TiCu 4 . Finally, a thick transition zone with non-defined borders can be seen between the fourth TiCu 4 + Ag s layer and the free solidified Ag-Cu(Ti) alloy. According to SEM, OM and EPMA, three phases were present at 790 °C in this transition zone: nearly eutectic Ag-Cu(Ti) alloy, a solid silver-rich phase and globular TiCu 4 crystals. The reaction layer sequence taking all these features into account can be expressed as: αTi / Ti 2 Cu / TiCu / Ti 3 Cu 4 + Ag s / TiCu 4 + Ag s / TiCu 4 + Ag s + L / L, where L designates a Ag s + Cu s nearly eutectic mixture resulting from the solidification of the liquid Ag-Cu(Ti) alloy (Cu s and TiCu 4 phases cannot be distinguished by SEM in Fig. 1 but are easily identifiable by their colour when observed by OM). It is to note that in the Ag-Cu-Ti system, no ternary eutectic transformation exists involving the four phases Ag s , Cu s , TiCu 4 and L. Instead, the liquid L reacts with TiCu 4 below 783 °C (transitory ternary peritectic reaction) to give solid copper and solid silver. This is the reason why although the liquid L is saturated in Ti, no TiCu 4 crystal precipitates up to complete solidification at 780 °C.

2. Ag-40 at.% Cu liquid without Ti

When no titanium was initially present in the Ag-Cu alloy, the major phenomenon seen was the dissolution of that element in the liquid alloy. After 50 min of heating at 790 °C, a depression with a depth attaining 60 µm was observable at the surface of the titanium plate, where a direct contact with the liquid alloy was established. Taking into account the area of the exposed surface (2 mm x 2 mm) and the volume of the liquid (∼80 mm 3 ), the increase of the Ti content in the drop was estimated to be about 0.2 at.%, which is only a tenth of the saturation value at 790 °C. Regardless of the heat-treatment time (from 2 min to 50 min), the reaction zone formed under these conditions at the liquid/solid interface remained very thin: its total thickness was in the range 3-6 µm. Fig. 2 shows that this zone consists of four different layers. Each of them is dense, continuous and exhibits a single-phased appearance, apart from some rare places in the outer layer where silver-rich inclusions are observable. The counting rates obtained for Ag, Cu and Ti by EPMA in the thickest parts of each layer were characteristic for single-phased compounds arranged according to the reaction layer sequence:

αTi / Ti 2 Cu / TiCu / Ti 3 Cu 4 / TiCu 4 / L', where L' designates a Ag s + Cu s mixture resulting from the solidification of the liquid alloy.

Fig. 2. Thin reaction zone formed between solid Ti and titanium-free Ag-40 at.% Cu liquid alloy after 50 min of heating at 790 °C

Reaction layer sequences and phase equilibria

Fig. 3 shows the isothermal section of the Ag-Cu-Ti phase diagram at 790 °C, calculated with the Thermo-calc software by using the already developed thermodynamic descriptions of the Ag-Cu [START_REF] Kattner | Phase Diagrams for Lead-free Solder Alloys[END_REF], Ag-Ti [START_REF] Li | Experimental Study and Thermodynamic Assessment of the Ag-Ti System[END_REF] and Cu-Ti [START_REF] Kumar | Thermodynamic Optimisation of the Cu-Ti System[END_REF] binaries. Extension to the third-order Ag-Cu-Ti system was made by using the assessment by Arroyave [START_REF] Arroyave | Thermodynamics and Kinetics of Ceramic/Metal Interfacial Interactions[END_REF]. In this section, which is in good agreement with previous literature data [START_REF] Eremenko | The liquidus Surface of the System Ti-Cu-Ag[END_REF][START_REF] Kubaschewski | Ternary Alloys[END_REF][START_REF] Villars | Handbook of Ternary Alloy Phase Diagrams[END_REF][START_REF]Ag-Cu-Ti (Silver-Copper-Titanium), Non-Ferrous Metal Ternary Systems. Selected Soldering and Brazing Systems: Phase Diagrams[END_REF], two reaction paths have been drawn to represent each of the two layer sequences. Both reaction paths comply with the mass balance, provided that the liquid is slightly enriched in Ag near the liquid/solid interface. For both sequences, it has moreover been verified by using the Thermo-Calc software that the activities of Ti and Cu regularly decrease and increase, respectively, when moving through the reaction zone from Ti to the liquid alloy as shown in Fig. 4. Accordingly, Ti and Cu atoms can intrinsically diffuse in opposite directions down their activity gradient [START_REF] Van Loo | Multiphase Diffusion in Binary and Ternary Solid-State Systems[END_REF]. Finally, local phase equilibria are established at all interfaces [START_REF] Yokokawa | Generalized Chemical Potential Diagram and Its Applications to Chemical Reactions at Interfaces between Dissimilar Materials[END_REF], except at the Ti 3 Cu 4 /TiCu 4 interface where the stable compound Ti 2 Cu 3 should be present. In fact, it has been shown in previous work that in the solid state, Ti 2 Cu 3 grows at a much slower rate than TiCu 4 and Ti 3 Cu 4 [START_REF] Andrieux | Details on the Formation of Ti 2 Cu 3 in the Ag-Cu-Ti System in the Temperature Range 790 to 860[END_REF].

The two reaction paths drawn in Fig. 3 can thus be considered as diffusion paths, it being understood that due to a very slow rate of growth, the compound Ti 2 Cu 3 is transiently "missing".

Discussion

Formation of the Ti/Ag-Cu reaction zones of the second type (Fig. 2) can be described as follows. As soon as the eutectic Ag-Cu alloy melts (M.P. 780 °C), the Ti base metal begins to dissolve in the liquid. When the Ti activity in the liquid is sufficiently high near the Ti surface, the compound TiCu 4 nucleates on Ti and rapidly forms a continuous reaction layer of small crystals that isolates Ti from the liquid as previously reported [START_REF] Andrieux | Details on the Formation of Ti 2 Cu 3 in the Ag-Cu-Ti System in the Temperature Range 790 to 860[END_REF]. Other Ti-Cu compounds can then nucleate and grow by solid state reaction-diffusion between Ti and the TiCu 4 layer. A reaction zone is thus formed (Fig. 4a), consisting of four dense and singlephased layers of compounds arranged according to the diffusion path B in the Ag-Cu-Ti phase diagram (Fig. 3). Obviously, the four layers tend to grow by solid state volume diffusion. In such an instance where growth by solid state diffusion and dissolution are competing, the reaction zone rapidly attains a constant thickness, as described in detail by Dibkov [START_REF] Dibkov | Growth Kinetics of Chemical Compound Layers[END_REF]. This explains why, although growing by solid state diffusion, the reaction zone does not increase in thickness. As a consequence, the Ti and Cu chemical activity gradients through the reaction zone remain steep (Fig. 4a), which might be the key for keeping up nearly planar interfaces. When the Ag-Cu liquid is fully saturated in titanium, dissolution of this element does not proceed and the different layers of the interface reaction zone can grow according to two main processes: solid state diffusion of Ti and Cu atoms for all layers and precipitation from a liquid locally over-saturated in Ti for the TiCu 4 crystals of the outer layer. As TiCu 4 precipitation proceeds at a fast rate, the liquid in contact with the outer layer is depleted in copper. This gives rise to the formation of some solid silver grains. Then, new TiCu 4 crystals precipitate on these grains from the liquid and new silver grains appear. Finally, a thick outer layer consisting of TiCu 4 crystals embedded in a solid silver-based matrix is formed after a few minutes of reaction at 790 °C. At this point, further precipitation of TiCu 4 and of silverrich grains implies the solid state diffusion of titanium through the growing reaction layers. Accordingly, as the thickness of the reaction layers increases, the rates of precipitation of TiCu 4 and Ag decrease. This would explain why the size and shape of the TiCu 4 crystals (and also of the silver-rich grains) gradually change from small and elongated to larger and more globular. As the thickness of the reaction zone continues to increase, another process, implying long-range crossed diffusion of Cu and Ti in the solid state, takes place. This process is the chemical conversion of the oldest TiCu 4 particles present in the inner part of the reaction zone into Ti 3 Cu 4 . Development of these different processes gives the liquid-solid reaction zone the morphological features schemed in Fig. 4b. It should be noted that, because of the saturation of the liquid in Ti and also of the large thickness of the two outer layers, the chemical activity gradients for Ti and Cu at the solid-liquid interface tend to be flat (Fig. 4b). This might be a reason why the solid/liquid boundary is not well-defined.

Let us now consider the formation of joints when Ti plates are brazed at 780-800 °C with Ag-40 at.% Cu alloy. For that purpose, 50-100 µm thick Ag-Cu ribbons are generally used, which means that solid titanium in excess reacts with a limited amount of Ag-Cu. The liquid alloy is thus rapidly saturated in Ti and the Ti/alloy reaction zone rapidly changes from type 2 to type 1. Fig. 5 shows a Ti/Ag-40 at.% Cu /Ti joint brazed at 790 °C for a short period of time (3 min) with an Ag-Cu ribbon. The particular morphology of this joint is the result of a time sequence of four of the above-described processes, namely Ti dissolution, growth by solid state interdiffusion, fast rate precipitation of TiCu 4 , and isothermal solidification of Ag. It is worth emphasizing that the four processes develop at a fast rate. Indeed, isothermal solidification of the braze alloy is completed after only 3 min of brazing. This implies fast rate migration of titanium and copper through solid silver. For lattice diffusion of Ti in solid Ag at 790 °C, a diffusivity of the order of 1.5 10 -9 cm 2 s -1 is reported [START_REF] Neumann | Application of the Modified Electrostatic Model to Diffusion of Transition Metals in Noble Metals[END_REF][START_REF] Zou | Formation of Titanium Nitride by Annealing Ag/Ti Structures in Ammonia Ambient[END_REF]. This value is too low to account for an isothermal solidification of the joint in less than 3 min. On the other hand, a diffusivity of the order of 5 10 -6 cm 2 s -1 is reported for grain boundary diffusion of Ti in Ag at 790 °C [START_REF] Zou | Formation of Titanium Nitride by Annealing Ag/Ti Structures in Ammonia Ambient[END_REF]. Such a high grain boundary diffusivity is consistent with the fast rate isothermal solidification observed, and this all the more than silver grains tend to grow in the form of joining columns (Fig. 5). The fact that a silver-rich zone without TiCu 4 crystals exists in the central part of the joint would indicate that Cu diffuses slightly faster than Ti. Upon further annealing (up to 15 hours) at 790 °C, the crystals of TiCu 4 embedded in the Ag-based matrix will be slowly converted by solid state reaction-diffusion first into Ti 3 Cu 4 and then into TiCu. This will increase the re-melting temperature of the brazed joint step by step from 780 to about 910 °C [START_REF] Dezellus | Transient Liquid Phase Bonding of Titanium to Aluminium Nitride[END_REF].

Conclusion

Interface reaction zones of two different types are formed by chemical interaction at 780-800 °C between solid titanium and liquid Ag-40 at.% Cu eutectic alloy: -as long as titanium is in contact with a pure Ag-Cu alloy, a thin reaction zone (3-6 µm in thickness) is present at the liquid/ solid interface. This zone consists of four continuous layers arranged according to the sequence αTi / Ti 2 Cu / TiCu / Ti 3 Cu 4 / TiCu 4 / L. The layer sequence is consistent with a solid state growth process, with Ti and Cu diffusing down their activity gradients. The fact that the reaction zone thickness does not change when the reaction time increases is explained by the fast dissolution rate of titanium in the liquid alloy; -when the Ag-Cu alloy is saturated in titanium (about 1.88 at.% Ti at 790 °C), dissolution no longer proceeds and a reaction zone is formed with a thickness increasing with the reaction time. In that case, the reaction layer sequence is: αTi / Ti 2 Cu / TiCu / Ti 3 Cu 4 + Ag s / TiCu 4 + Ag s / TiCu 4 + Ag s + L / L. In addition to Ti dissolution and Ti-Cu solid state interdiffusion, three other elementary reaction processes are identified. These are the nucleation and growth from the liquid of TiCu 4 particles, the isothermal solidification of silver, and the chemical conversion of Cu-Ti compounds by reaction-diffusion in the solid state. A time sequence of these five elementary reaction processes is finally proposed to account for the morphology and composition of Ti/ Ag-Cu/ Ti joints brazed at 780-800 °C.

Fig. 1 .

 1 Fig. 1. Thick reaction zone formed between solid Ti and liquid Ag-40 at.% Cu alloy saturated in Ti (1.88 at.%Ti) after 50 min of heating at 790 °C: (a) general view; (b) zoom near the Ti surface

Fig. 3 .

 3 Fig. 3. The Ag-Cu-Ti ternary phase diagram as calculated at 790 °C. Two diffusion paths have been drawn corresponding to the reaction zones formed in Ag-40 at.% Cu alloys saturated in Ti (path A) and without Ti (path B)

Fig. 4 .

 4 Fig. 4. Chemical activities of fcc-Cu and cph-Ti in the intermetallic compound layers formed by reaction at 790 °C of solid titanium with (a) pure and (b) Ti-saturated Ag-40 at.% Cu liquid alloy

Fig. 5 .

 5 Fig. 5. Ti/Ag-40 at.% Cu / Ti sandwich joint after 3 min of brazing at 790 °C

Table 1 .

 1 Ag, Cu and Ti contents found by EPMA in the reaction layers formed after 50 min at 790 °C between solid titanium and Ti-saturated Ag-40at.% Cu liquid alloy

			Ag	Cu	Ti	Phase
			(at.%)	(at.%)	(at.%)	attribution
	1 st layer		0 -2	34.0 -32.0	66.0	Ti 2 Cu
	2 nd layer		2.0 -3.3	48.0 -46.7	50.0	TiCu
	3 rd layer	dark phase bright phase	1.3 -1.9 ∼ 86	55.2-54.6 ∼ 9	43.5 ∼ 5	Ti 3 Cu 4 Ag
	4 th layer	dark phase bright phase	1.5 -2.0 ∼ 86	75-74.5 ∼ 14	23.5 -	TiCu 4 Ag


Julie Delon

Julien Salomon

Andrei Sobolevskii

AND Andrei Sobolevski

FAST TRANSPORT OPTIMIZATION FOR MONGE COSTS ON THE CIRCLE *

Keywords: Monge-Kantorovich problem, Monge cost, Aubry-Mather (weak KAM) theory AMS subject classifications. Primary 90C08; Secondary 68Q25, 90C25

published or not. The documents may come

Fast transport optimization for Monge costs on the circle

1. Introduction. The transport optimization problem, introduced by G. Monge in 1781 and shown by L. Kantorovich in 1942 to be an instance of linear programming, is a convex optimization problem with strong geometric features. A typical example is minimization of mean-square displacement between two given finite marginal measures supported on convex compacta in Euclidean space: in this case a solution is defined by gradient of a convex function that satisfies a suitable Monge-Ampère equation. Various generalizations of this result and rich bibliographies can be found, e.g., in the survey [START_REF] Gangbo | The geometry of optimal transportation[END_REF] or the recent monograph [START_REF] Villani | Optimal transport: Old and new[END_REF].

Further constraints on the two marginals or their supports may furnish the problem with useful additional convex structure. One way of making this statement quantitative is to consider the algorithmic complexity of the corresponding numerical transport optimization schemes. In particular when the two measures live on segments of straight lines, the optimal map is monotone and may be found by sorting, which takes O(n log n) operations when the data come in the form of discrete n-point histograms. If the input data are already sorted, this count falls to O(n).

The optimal transport is well understood also when the marginals live on a compact Riemannian manifold [START_REF] Feldman | Monge's transport problem on a Riemannian manifold[END_REF]; the existence and characterization of optimal map in the case of a flat torus and quadratic cost have been established a decade ago [START_REF] Cordero-Erausquin | Sur le transport de mesures périodiques[END_REF]. However, the algorithmics of even the simplest setting of the unit circle is no longer trivial, because the support of the measures is now oriented rather than ordered. A naive approach would require solving the problem for each of n different alignments of two n-point histograms, thereby involving O(n 2 ) operations.

In this paper we present an efficient algorithm of transport optimization on the circle, which is based on a novel analogy with the Aubry-Mather (weak KAM) theory in Lagrangian dynamics (see, e.g., [3,[START_REF] Fathi | Weak KAM Theorem in Lagrangian Dynamics[END_REF][START_REF] Knill | selected chapters in the calculus of variations[END_REF]). The key step is to lift the transport problem to the universal cover of the unit circle, rendering the marginals periodic and the cost of transport infinite. However, it still makes sense to look for those transport maps whose cost cannot be decreased by any local modification. Different locally optimal maps, which cannot be deformed into each other by any local rearrangement, form a family parameterized with an analogue of the rotation number in the Aubry-Mather theory. One can introduce a counterpart of the average Lagrangian, or Mather's α function in the Aubry-Mather theory, which turns out to be efficiently computable. As we show below, its minimization provides an efficient algorithm of transport optimization on the circle. The class of cost functions for which this theory works includes all costs with the Monge property, such as the quadratic cost or costs generated by natural Lagrangians with time-periodic potentials [START_REF] Bernard | Optimal mass transportation and Mather theory[END_REF][START_REF] Knill | selected chapters in the calculus of variations[END_REF].

Note that the problem of optimally matching circular distributions appears in a variety of applications. Important examples are provided by image processing and computer vision: image matching techniques for retrieval, classification, or stitching purposes [START_REF] Zhang | Local features and kernels for classification of texture and object categories: A comprehensive study[END_REF][START_REF] Brown | Multi-image matching using multi-scale oriented patches[END_REF] are often based on matching or clustering "descriptors" of local features [START_REF] Lowe | Distinctive image features from scale-invariant keypoints[END_REF], which typically consist of one or multiple histograms of gradient orientation. Similar issues arise in object pose estimation and pattern recognition [START_REF] Lowe | Distinctive image features from scale-invariant keypoints[END_REF][START_REF] Gangbo | Shape recognition via Wasserstein distance[END_REF]. Circular distributions also appear in a quite different context of analysis of color images, where hue is parameterized by polar angle. In all these applications, matching techniques must be robust to data quantization and noise and computationally effective, which is especially important with modern large image collections. These requirements are satisfied by the optimal value of a transport cost for a suitable cost function.

This paper is organized as follows. In §2 we give a specific but nontechnical overview of our results. After the basic definitions are given in §3, including that of locally optimal transport plans, in §4 we give an explicit description of the family of locally optimal transport plans: they are conjugate, in measure theoretic sense, to rotations of the unit circle (or equivalently to shifts of its universal cover). This result is in direct analogy with conjugacy to rotations in the one-dimensional Aubry-Mather theory [3]. As shown in §5, the average cost of a locally optimal transport plan is a convex function of the shift parameter. Moreover, the values of this function and its derivative are efficiently computable when the marginal measures are discrete, which enables us to present in §6 a fast algorithm for transport optimization on the circle. The same section contains results of a few numerical experiments. Finally a review of related work in the computer science literature is given in §7.

2. Informal overview. For two probability measures μ0 , μ1 on the unit circle T = R/Z and a given cost ĉ(x, ŷ) of transporting a unit mass from x to ŷ in T, the transport cost is defined as the inf of the quantity

Î(γ) = T×T ĉ(x, ŷ) γ(dx × dŷ). (2.1)
over the set of all couplings γ of the probability measures μ0 , μ1 (i.e., all measures on T × T with marginals μ0 , μ1 ). These couplings are usually called transport plans. Suppose that the cost function ĉ(•, •) on T × T is determined via the relation ĉ(x, ŷ) = inf c(x, y) by a function c(•, •) on R×R satisfying the condition c(x+1, y+1) = c(x, y) for all x, y; here inf is taken over all x, y whose projections to the unit circle coincide with x, ŷ. We lift the measures μ0 and μ1 to R, obtaining periodic locally

v u O F 1 F θ 1 -θ (F θ 1 ) -1 (v) F -1 0 (v) v Fig. 2.1.
Construction of the locally optimal transport plan γ θ .

finite measures µ 0 , µ 1 , and redefine γ to be their coupling on R × R. It is then convenient to replace the problem of minimizing the integral (2.1) with "minimization" of an integral

I(γ) = R×R c(x, y) γ(dx × dy). (2.2)
Although the latter integral is infinite, it still makes sense to look for transport plans γ minimizing I with respect to local modifications, i.e., to require that for any compactly supported signed measure δ of zero mass and finite total variation, the difference I(γ + δ) -I(γ), which is defined by a finite integral, be nonnegative. These locally optimal transport plans are the main object of this paper. Assume that the cost function c(x, y) satisfies the Monge condition (alternatively known as the continuous Monge property, see [1,[START_REF] Burkard | Perspectives of Monge properties in optimization[END_REF]):

c(x 1 , y 1 ) + c(x 2 , y 2 ) < c(x 1 , y 2 ) + c(x 2 , y 1 ) (2.3)
for all x 1 < x 2 and y 1 < y 2 . An example of such a cost function is |x -y| λ , where λ > 1; in this case the quantity MK λ (μ 0 , μ1 ) = (inf γ Î(γ)) 1/λ turns out to be a metric on the set of measures on the circle, referred to as the Monge-Kantorovich distance of order λ. The value λ = 1 can still be treated in the same framework as the limiting case λ → 1; it is sometimes called the Kantorovich-Rubinshtein metric or, in image processing literature, the Earth Mover's distance [START_REF] Rubner | The Earth Mover's Distance as a Metric for Image Retrieval[END_REF]. The Monge condition (2.3) implies that whenever under a transport plan the mutual order of any two elements of mass is reversed, the transport cost can be strictly reduced by exchanging their destinations. It follows that a locally minimal transport plan moves elements of mass monotonically, preserving their spatial order.

The whole set of locally optimal transport plans for a given pair of marginals µ 0 , µ 1 can be conveniently described using a construction represented in fig. 2.1. Let F 0 , F 1 be cumulative distribution functions of the measures µ 0 , µ 1 normalized so that F 0 (0) = F 1 (0) = 0. We shall regard graphs of F 0 , F 1 as continuous curves including, where necessary, the vertical segments corresponding to jumps of these functions (which are caused by atoms of µ 0 , µ 1 ). Each of these curves specifies a correspondence, F -1

0 or F -1
1 , between points of the vertical axis Ov, representing elements of mass, and points of the horizontal axis Ou, representing spatial locations, and maps the Lebesgue measure on Ov into µ 0 or µ 1 on Ou. This correspondence is monotone and defined everywhere except on an (at most countable) set of v values that correspond to vacua of the measure in the Ou axis.

Define now F θ 1 (u) = F 1 (u)-θ. Then (F θ 1 ) -1 represents a shift of the Ov axis by θ followed by an application of the correspondence F -1

1 , and still induces on the Ou axis the same measure µ 1 as F 1 . A transport plan γ θ that takes an element of mass represented by v from F -1 0 (v) to (F θ 1 ) -1 (v) is, by construction, a monotone coupling of µ 0 and µ 1 , and thus a locally optimal transport plan. Moreover, it is shown in §4 that all locally optimal transport plans can be obtained using this construction as the parameter θ runs over R.

Finally define the average cost C [F0,F1] (θ) of the plan γ θ per unit period:

C [F0,F1] (θ) = 1 0 c(F -1 0 (v), (F θ 1 ) -1 (v)) dv.
It is shown in §5 that the Monge condition implies convexity of C [F0,F1] (θ) and that the global minimum of this function in θ coincides with the minimum value of the transport cost on the unit circle (2.1).

When the marginals µ 0 , µ 1 are purely atomic with finite numbers n 0 and n 1 of atoms in each period, the function C becomes piecewise affine. In §6 we present an algorithm to approximate its minimum value to accuracy ǫ, using a binary search that takes O((n 0 + n 1 ) log(1/ǫ)) operations in the real number computing model. When masses of all atoms are rational numbers with the least common denominator M , this search returns an exact solution provided that ǫ < 1/M . This gives an O((n 0 + n 1 ) log M ) exact transport optimization algorithm on the circle.

Preliminaries.

Let T = R/Z be the unit circle, i.e., the segment [0, 1] with identified endpoints. By π : R → T denote the projection that takes points of the universal cover R to points of T.

The cost function.

A cost function is a real-valued function c(•, •) defined on the universal cover R of the circle T. We assume that it satisfies the Monge condition: for any x 1 < x 2 and y 1 < y 2 , c(x 1 , y 1 ) + c(x 2 , y 2 ) -c(x 1 , y 2 ) -c(x 2 , y 1 ) < 0.

(3.1)

Additionally c is assumed to be lower semicontinuous, to be invariant with respect to integer shifts, i.e.,

c(x + 1, y + 1) = c(x, y) (3.2)
for all x, y, and to grow uniformly as |x -y| → ∞: for any P there exists a finite R(P ) ≥ 0 such that c(x, y) ≥ P whenever |x -y| ≥ R(P ). (3.3) Note that the latter condition implies that the lower semicontinuous function c is bounded from below (and guarantees that the minima in a number of formulas below are attained).

Note that the Monge condition (3.1) holds for any twice continuously differentiable function c such that ∂ 2 c(x, y)/∂x ∂y < 0. If the cost function depends only on x -y, this reduces to a convexity condition: -∂ 2 c(x -y)/∂x ∂y = c ′′ (x -y) > 0. In particular, all the above conditions are satisfied for the function c(x, y) = |x -y| λ , which appears in the definition of the Monge-Kantorovich distance with λ > 1, and, more generally, for any function of the form c(x-y)+f (x)+g(y) with strictly convex c and periodic f and g.

For a cost function c satisfying all the above conditions, the cost of transporting a unit mass from x to ŷ on the circle is defined as ĉ(x, ŷ) = inf c(x, y), where x, ŷ are points of T and inf is taken over all x, y in R such that πx = x and πy = ŷ. Using the integer shift invariance, this definition can be lifted to the universal cover as ĉ(x, y) = inf k∈Z c(x, y + k).

Condition (3.1) is all that is needed in §4, which is concerned with locally optimal transport plans on R. Conditions (3.2), (3.3) come into play in §5, which deals with transport optimization on the circle.

Distribution functions.

For a given locally finite measure µ on R define its distribution function F µ by

F µ (0) = 0, F µ (x) = µ((0, x]) for x > 0, F µ (x) = -µ((x, 0]) for x < 0. (3.4) Then µ((x 1 , x 2 ]) = F µ (x 2 ) -F µ (x 1 ) whenever x 1 < x 2
, and this identity also holds for any function that differs from F µ by an additive constant (the normalization F µ (0) = 0 is arbitrary). When µ is periodic with unit mass in each period, it follows that for all

x in R F µ (x + 1) = F µ (x) + 1. (3.5)
The inverse of a distribution function F µ is defined by are right-continuous. Discontinuities of F µ correspond to atoms of µ and discontinuities of its inverse, to "vacua" of µ, i.e., to intervals of zero µ measure.

F -1 µ (y) = inf{x : y < F µ (x)} = sup{x : y ≥ F µ (x)}. ( 3
For a distribution function F µ define its complete graph to be the continuous curve formed by the union of the graph of F µ and the vertical segments corresponding to jumps of F µ . Accordingly, by a slight abuse of notation let F µ ({x}) denote the set [F µ (x -0), F µ (x)] (warning: F µ ({x}) ⊇ {F µ (x)}) and let F µ (A) = x∈A F µ ({x}) for any set A.

Local properties of transport plans.

Let μ0 , μ1 be two finite positive measures of unit total mass on T and µ 0 , µ 1 their liftings to the universal cover R, i.e., periodic measures such that µ i (A) = μi (πA), i = 0, 1, for any Borel set A that fits inside one period. Periodicity of measures here means that µ(A + n) = µ(A) for any integer n and any Borel A, where A + n = {x + n : x ∈ A}.

Definition 3.1. A (locally finite)1 transport plan with marginals µ 0 and µ 1 is a locally finite measure γ on R × R such that (i) for any x in R the supports of measures γ((-∞, x] × •) and γ(• × (-∞, x]) are bounded from above and the supports of measures γ

((x, ∞) × •), γ(• × (x, ∞)) are bounded from below; (ii) γ(A × R) = µ 0 (A) and γ(R × B) = µ 1 (B) for any Borel sets A, B.
The quantity γ(A × B) is the amount of mass transferred from A to B under the transport plan γ. Condition (1) implies that the mass supported on any bounded interval gets redistributed over a bounded set (indeed, a bounded interval is the intersection of two half-lines), but is somewhat stronger. Definition 3.2. A local modification of the locally finite transport plan γ is a transport plan γ ′ such that γ and γ ′ have the same marginals and γ ′ -γ is a compactly supported finite signed measure. A local modification is called cost-reducing if

c(x, y) (γ ′ (dx × dy) -γ(dx × dy)) < 0.
A locally finite transport plan γ is said to be locally optimal with respect to the cost function c or c-locally optimal if it has no cost-reducing local modifications.

4. Conjugate transport plans and shifts. Let U 0 , U 1 be two copies of R equipped with positive periodic measures µ 0 , µ 1 whose distribution functions F 0 , F 1 satisfy (3.5), so that all intervals of unit length have unit mass. Let furthermore V 0 , V 1 be two other copies of R equipped with the uniform (Lebesgue) measure.

4.1. Normal plans and conjugation. We introduce the following terminology:

Definition 4.1. A locally finite transport plan ν on V 0 × V 1 with uniform marginals is called normal. Definition 4.2. For a normal transport plan ν its conjugate transport plan ν [F0,F1] is a transport plan on U 0 × U 1 such that for any Borel sets A, B ν [F0,F1] (A × B) = ν(F 0 (A) × F 1 (B)). ( 4.1)
Lemma 4.3. For a normal transport plan ν its conjugate ν [F0,F1] is a locally finite transport plan on U 0 × U 1 with marginals µ 0 , µ 1 .

Proof. Since distribution functions F 0 , F 1 and their inverses preserve boundedness, condition (1) of Definition 3.1 is fulfilled. Definition 4.2, condition (2) of Definition 3.1, and formula (3.4) together imply that F1] satisfies condition (2) of Definition 3.1 on intervals and therefore on all Borel sets. Lemma 4.4. For any transport plan γ on U 0 × U 1 with marginals µ 0 and µ 1 there exists a normal transport plan ν such that γ is conjugate to ν: γ = ν [F0,F1] .

ν [F0,F1] ((u 1 , u 2 ] × U 1 ) = ν(F 0 ((u 1 , u 2 ]) × F 1 (U 1 )) = ν([F 0 (u 1 ), F 0 (u 2 )] × V 1 ) = F 0 (u 2 ) -F 0 (u 1 ) = µ 0 ((u 1 , u 2 ]). Similarly ν [F0,F1] (U 0 × (u 1 , u 2 ]) = µ 1 ((u 1 , u 2 ]). Thus ν [F0,
Proof. For non-atomic measures µ 0 and µ 1 the required transport plan is given by the formula

ν(A × B) = γ(F -1 0 (A) × F -1 1 (B))
, which is dual to (4.1). However if, e.g., µ 0 has an atom, then the function F -1 0 is constant over a certain interval and maps any subset A of this interval into one point of fixed positive measure in U 0 , so information on the true Lebesgue measure of A is lost. In this case extra care has to be taken.

Recall that a locally finite measure has at most a countable set of atoms. Let atoms of µ 0 be located in [START_REF] Sobolevski | then θ is the required minimum; stop[END_REF]1]

at points u 1 , u 2 , . . . with masses m 1 , m 2 , . . . . Since γ({u i } × U 1 ) = µ 0 ({u i }) = m i > 0, there exists a conditional probability measure ρ(• | u i ) = γ({u i } × •)/m i . For a set A ⊂ (0, 1] define a "residue" transport plan γ(A × B) = γ(A × B) -i m i δ ui (A) ρ(B | u i ),
where δ u is the Dirac unit mass measure on U 0 concentrated at u, and extend γ to general A using periodicity. We thus remove from γ the part of γ whose projection to the first factor is atomic. Define a transport plan κ on V 0 × U 1 by

κ(C × B) = i λ(C ∩ F 0 ({u i })) ρ(B | u i ) + γ(F -1 0 (C) × B),
where C is a Borel set in V 0 and λ(•) denotes the Lebesgue measure in V 0 . Clearly κ(F 0 (A)×B) = γ(A×B). Repeating this construction for the second factor, with κ in place of γ, we get a normal transport plan ν such that γ(

A × B) = ν(F 0 (A) × F 1 (B)).
Since we are ultimately interested in transport optimization with marginals µ 0 , µ 1 rather than with uniform marginals, two normal transport plans ν 1 , ν 2 will be called equivalent if they have the same conjugate. Two different normal transport plans can only be equivalent if one or both measures µ 0 or µ 1 have atoms, causing loss of information on the structure of ν in segments corresponding to these atoms. The proof of Lemma 4.4 gives a specific representative of this equivalence class of normal plans.

4.2.

Locally optimal normal transport plans are shifts. Fix a cost function c : U 0 × U 1 → R that satisfies the Monge condition (3.1) and define

c [F0,F1] (v 0 , v 1 ) = c F -1 0 (v 0 ), F -1 1 (v 1 ) . (4.2)
For non-atomic measures µ 0 , µ 1 , it satisfies the Monge condition Here and below, variables u, u ′ , u 0 , u 1 , . . . are assumed to take values in U 0 or U 1 and variables v, v ′ , v 0 , v 1 , . . . , w, w ′ , . . . , in V 0 or V 1 .

c [F0,F1] (v ′ , w ′ ) + c [F0,F1] (v ′′ , w ′′ ) -c [F0,F1] (v ′ , w ′′ ) -c [F0,F1] (v ′′ , w ′ ) < 0 whenever v ′ < v ′′ and
Lemma 4.5. A transport plan γ on U 0 × U 1 with marginals µ 0 , µ 1 is c-locally optimal if and only if it is conjugate to a c [F0,F1] -locally optimal normal transport plan ν. In particular, all normal transport plans with the same locally optimal conjugate are locally optimal.

Proof. Note that ν ′ -ν is compactly supported if and only if the difference of the respective conjugates γ ′ -γ is compactly supported. The rest of the proof follows from the identity

c [F0,F1] (v 1 , v 2 ) ν ′ (dv 1 × dv 2 ) -ν(dv 1 × dv 2 ) = c(u 1 , u 2 ) γ ′ (du 1 × du 2 ) -γ(du 1 × du 2 )
established by the change of variables v 1 = F 0 (u 1 ), v 2 = F 1 (u 2 ) (here jumps of the distribution functions are harmless because c [F0,F1] is constant over respective ranges of its variables).

Transport optimization with marginals µ 0 , µ 1 is thus reduced to a conjugate problem involving uniform marginals and the cost c [F0,F1] . It turns out that any c [F0,F1] -optimal normal transport plan must be supported on a graph of a monotone function, and due to uniformity of marginals this function can only be a shift by a suitable real increment θ. More precisely, the following holds:

Theorem 4.6. Let µ 0 , µ 1 be two periodic positive measures defined respectively on U 0 , U 1 with unit mass in each period and let F i : U i → V i , i = 0, 1, be their distribution functions. Then any c [F0,F1] -locally optimal normal transport plan on V 0 × V 1 is equivalent to a normal transport plan ν θ with supp ν θ = {(v, w) : w = v + θ}, and conversely ν θ is c [F0,F1] -locally optimal for any real θ. All c-locally optimal transport plans on U 0 × U 1 with marginals µ 0 , µ 1 are of the form γ θ = (ν θ ) [F0,F1] .

The proof, divided into a series of lemmas, is based on the classical argument: a nonoptimal transport plan can be modified by "swapping" pieces of mass to render its support monotone while decreasing its cost. This argument, carried out for plans with uniform marginals on V 0 × V 1 , is combined with the observation that a monotonicaly supported plan with uniform marginals can only be a shift. Then Lemma 4.5 is used to extend this result to transport plans on U 0 × U 1 .

Throughout the proof fix a normal transport plan ν and define on

V 0 × V 1 the functions r ν (v, w) = ν((-∞, v] × (w, ∞)), l ν (v, w) = ν((v, ∞) × (-∞, w]). (4.3)
To explain the notation r ν , l n u observe that, e.g., r ν (v, w) is the amount of mass that is located initially to the left of v and goes to the right of w. Lemma 4.7. The function r ν (resp. l ν ) is continuous and monotonically increasing in its first (second) argument and is continuous and monotonically decreasing in its second (first) argument, while the other argument is kept fixed.

Proof. Monotonicity is obvious from (4.3). To prove continuity observe that the second marginal of ν is uniform, which together with positivity of all involved measures implies that in the decomposition

ν(V 0 × • ) = ν((-∞, v] × • ) + ν((v, ∞) × • ),
both measures in the right-hand side cannot have atoms. This implies continuity of r ν , l ν with respect to the second argument. A similar proof holds for the first argument.

Lemma 4.8. For any v there exist w ν (v) and m ν (v) ≥ 0 such that

r ν (v, w ν (v)) = l ν (v, w ν (v)) = m ν (v). ( 4

.4)

The correspondence v → w ν (v) is monotone:

w ν (v 1 ) ≤ w ν (v 2 ) for v 1 < v 2 . Proof. Clearly r ν (v, -∞) = ∞, r ν (v, ∞) = 0, l ν (v, -∞) = 0, l ν (v, ∞) = ∞.
The continuity of the functions r ν (v, •), l ν (v, •) in the second argument for a fixed v implies that their graphs intersect at some point (w ν (v), m ν (v)), which satisfies (4.4). Should the equality r ν (v, w) = l ν (v, w) hold on a segment [w ′ , w ′′ ], we set w ν (v) to its left endpoint w ′ ; this situation, however, will be ruled out by the corollary to Lemma 4.10 below. Monotonicity of w ν (v) follows from monotonicity of r ν (•, w), l ν (•, w) in the first argument for a fixed w: indeed, for

v 2 > v 1 the equality r ν (v 2 , w) = l ν (v 2 , w) is impossible for w < w ν (v 1 ) because for such w we have r ν (v 2 , w) > r ν (v 1 , w ν (v)) = l ν (v 1 , w ν (v)) > l ν (v 2 , w).
Equalities (4.4) mean that the same amount of mass m ν (v) goes under the plan ν from the left of v to the right of w ν (v) and from the right of v to the left of w ν (v).

We are now in position to use the Monge condition and show that this amount can be reduced to zero by modifying the transport plan locally without a cost increase.

Lemma 4.9. For any v there exists a local modification ν v of ν such that w νv (v) = w ν (v) (with w ν defined as in Lemma 4.8), m νv (v) = 0, and ν v is either cost-reducing in the sense of Definition 3.2 or is equivalent to ν.

Proof. Let w = w ν (v) and m = m ν (v). If m = 0, there is nothing to prove. Suppose that m > 0 and define

w -= sup{w ′ : l ν (v, w ′ ) = 0}, w + = inf{w ′ : r ν (v, w ′ ) = 0}, v -= sup{v ′ : r ν (v ′ , w) = 0}, v + = inf{v ′ : l ν (v ′ , w) = 0}.
By local finiteness of the transport plan ν all these quantities are finite. Since m > 0, continuity of r ν , l ν implies that the inequalities w -< w < w + and v -< v < v + are strict. Consider the measures

ρ -(•) = ν( • × (w, w + )) on (v -, v), ρ + (•) = ν( • × (w -, w)) on (v, v + ), σ -(•) = ν((v, v + ) × • ) on (w -, w), σ + (•) = ν((v -, v) × • ) on (w, w + ).
Equalities (4.4) mean that all these measures have the same positive total mass m. Note that the Lebesgue measures of intervals (v -, v), (v, v + ), (w -, w), and (w, w + ) may be greater than m, because some mass in these intervals may come from or go to elsewhere. The functions r w (

•) = r ν (•, w), l v (•) = l ν (v, •) are monotonically increasing and r v (•) = r ν (v, •), l w (•) = l ν (•, w) are monotonically decreasing, with their inverses r -1 w , l -1 v , r -1 v , l -1
w defined everywhere except on an at most countable set of points. These functions may be regarded as a kind of distribution functions for the measures ρ -, σ -, σ + , ρ + respectively, mapping them to the Lebesgue measure on (0, m).

Under the plan ν, mass m is sent from (v -, v) to (w, w + ) and from (v, v + ) to (w -, w). We now construct a local modification ν v of the transport plan ν that moves mass m from the interval (v -, v) to (w -, w) and from (v, v + ) to (w, w + ), and show that it is cost-reducing unless measures µ 0 , µ 1 have atoms corresponding to the intervals under consideration.

Observe first that the normal plan ν induces two transport plans τ r , τ l that map measures ρ -to σ + and ρ + to σ -correspondingly:

τ r (A × B) = ν(A ∩ (-∞, v) × B ∩ (w, +∞)) = ν(A ∩ (v -, v) × B ∩ (w, w + )), τ l (A × B) = ν(A ∩ (v, +∞) × B ∩ (-∞, w)) = ν(A ∩ (v, v + ) × B ∩ (w -, w)),
where A ⊂ (v -, v + ), B ⊂ (w -, w + ) are two arbitrary Borel sets and the ∩ operation takes precedence over ×. By an argument similar to the proof of Lemma 4.4, there exist two transport plans χ r and χ l mapping the Lebesgue measure on (0, m) respectively to σ + , σ -and such that

τ r (A × B) = χ r r w (A ∩ (v -, v)) × B ∩ (w, w + ) , τ l (A × B) = χ l l w (A ∩ (v, v + )) × B ∩ (w -, w) .
Define now two transport plans τl , τr that send mass elements to the same destinations but from interchanged origins:

τr (A × B) = χ r l w (A ∩ (v, v + )) × B ∩ (w, w + ) , τl (A × B) = χ l r w (A ∩ (v -, v)) × B ∩ (w -, w) .

This enables us to define

ν v (A × B) = ν(A × B) -τ r (A × B) -τ l (A × B) + τr (A × B) + τl (A × B).
Since τ r (A×R) = τr (A×R) = ρ -(A) etc., the transport plan ν v has the same uniform marginals as ν, i.e., it is a local modification of ν. Observe furthermore that by the construction of ν v no mass is moved under this plan from the left-hand side of v to the right-hand side of w and inversely, i.e., that m νv (v) = 0.

It remains to show that ν v is either a cost-reducing modification of ν or equivalent to it. By the disintegration lemma (see, e.g., [START_REF] Ambrosio | Gradient Flows: In Metric Spaces and in the Space of Probability Measures[END_REF]) we can write χ r (dα

× dw ′ ) = dα dG r (w ′ | α) and χ l (dα × dw ′ ) = dα dG l (w ′ | α), where G r ( • | α) (resp. G l ( • | α))
are distribution functions of probability measures defined on [w, w + ] (resp. [w -, w]) for almost all 0 < α < m. Denote their respective inverses by

G -1 r (• | α), G -1 l (• | α) and observe that w -≤ G -1 l (β ′ | α) ≤ w ≤ G -1 r (β ′′ | α) ≤ w + for any β ′ , β ′′ . Thus c(v ′ , w ′ ) τ r (dv ′ × dw ′ ) = c(r -1 w (α), w ′ ) χ r (dα × dw ′ ) = m 0 dα c(r -1 w (α), w ′ ) dG r (w ′ | α) = m 0 dα 1 0 dβ c(r -1 w (α), G -1 r (β | α)),
where we write c instead of c [F0,F1] to lighten notation, and similarly

c(v ′ , w ′ ) τ l (dv ′ × dw ′ ) = m 0 dα 1 0 dβ c(l -1 w (α), G -1 l (β | α)), c(v ′ , w ′ ) τr (dv ′ × dw ′ ) = m 0 dα 1 0 dβ c(l -1 w (α), G -1 r (β | α)), c(v ′ , w ′ ) τl (dv ′ × dw ′ ) = m 0 dα 1 0 dβ c(r -1 w (α), G -1 l (β | α)).
The integral in Definition 3.2 now takes the form

c(v ′ , w ′ ) (ν v (dv ′ × dw ′ ) -ν(dv ′ × dw ′ )) = c(v ′ , w ′ ) -τ r (dv ′ × dw ′ ) -τ l (dv ′ × dw ′ ) + τr (dv ′ × dw ′ ) + τl (dv ′ × dw ′ ) = m 0 dα 1 0 dβ -c(r -1 w (α), G -1 r (β | α)) -c(l -1 w (α), G -1 l (β | α)) + c(l -1 w (α), G -1 r (β | α)) + c(r -1 w (α), G -1 l (β | α)) . As r -1 w (α) ≤ v ≤ l -1 w (α) and G -1 l (β | α) ≤ w ≤ G -1 r (β | α)
for all α, β, the Monge condition (3.1) implies that either the value of this integral is negative or the function c (i.e., c [F0,F1] ) is constant in at least one of its arguments. In the former case the transport plan ν v is a cost-reducing local modification of ν; in the latter case ν v is equivalent to ν.

Lemma 4.10. For any v ′ < v ′′ there exists a local modification

ν v ′ ,v ′′ of ν such that w ν v ′ ,v ′′ (v) = w ν (v) for v ′ ≤ v ≤ v ′′ , m ν v ′ ,v ′′ (v ′ ) = m ν v ′ ,v ′′ (v ′′ ) = 0, and in the strip v ′ ≤ v ≤ v ′′ the support of ν v ′ ,v ′′ coincides

with the complete graph of the monotone function w ν (•).

Proof. Let {v i } be a dense countable subset of [v ′ , v ′′ ] including its endpoints. Set ν 0 = ν and define ν i recursively to be the local modification of ν i-1 given by the previous lemma and such that w νi (v i ) = w ν (v i ) and m νi (v i ) = 0. Then all ν i are either cost-reducing or equivalent to ν and w νj (v i ) = w ν (v i ), m νj (v i ) = 0 for all j > i. Indeed, denote w i = w νi (v i ) and observe that if e.g. v j > v i , then, as m νi (v i ) = r νi (v i , w i ) = 0, mass from (-∞, v i ] does not appear to the right of w i and so does not contribute to the balance of mass around w j . Therefore for any j the possible modification of ν j-1 is local to the interval (v i ′ , v i ′′ ), where v i ′ = max{v i : i < j, v i < v j } and v i ′′ = min{v i : i < j, v i > v j } (with max and min of empty set defined as v ′ and v ′′ ). Thus there is a well-defined limit normal transport plan ν ∞ that is either a cost-reducing local modification or equivalent to ν and is such that, by continuity of the functions r ν∞ and l ν∞ in the first argument, m ν∞ (v) vanishes everywhere on [v ′ , v ′′ ].

Consider now the function w ν∞ (•), which coincides with w ν (•) on a dense subset of [v ′ , v ′′ ], so that their complete graphs coincide. For any quadrant of the form (-∞, v 0 ) × (w 0 , ∞) such that w 0 > w ν∞ (v 0 ), monotonicity of r ν∞ in the second argument implies that

0 ≤ ν ∞ ((-∞, v 0 ) × (w 0 , ∞)) = r ν∞ (v 0 , w 0 ) ≤ r ν∞ (v 0 , w ν∞ (v 0 )) = 0, i.e., ν ∞ ((-∞, v 0 ) × (w 0 , ∞)) = 0. Similarly ν ∞ ((v 0 , ∞) × (-∞, w 0 )) = 0
for any quadrant with w 0 < w ν∞ (v 0 ). The union of all such quadrants is the complement of the complete graph of the function v → w ν (v); this implies that ν ∞ is supported thereon.

Corollary 4.11. For any normal transport plan ν there exists a real number

θ ν such that w ν (v) = v + θ ν . Proof. It is enough to show that w ν (v ′ ) -v ′ = w ν (v ′′ ) -v ′′ for all v ′ , v ′′ . Let v ′ < v ′′
and ν v ′ ,v ′′ be the local modification constructed in the previous lemma.

Since it has uniform marginals and monotone support, we have

w ν (v ′′ ) -w ν (v ′ ) = ν v ′ ,v ′′ ((v ′ , v ′′ ) × (w ν (v ′ ), w ν (v ′′ ))) = v ′′ -v ′ , which completes the proof.
We call the parameter θ ν the rotation number of the normal transport plan ν. Definition 4.12. A normal transport plan consisting of a uniform measure supported on the line {(v, w) : w = v + θ} is called a shift and denoted by ν θ . Lemma 4.13. For any θ the shift ν θ is c [F0,F1] -locally optimal. Proof. Let ν be a local modification of ν θ such that the signed measure ν θν is supported in (v ′ , v ′′ ) × (w ′ , w ′′ ). Let νv ′ ,v ′′ be a local modification of ν constructed in Lemma 4.10; it coincides with ν θ over v ′ < v < v ′′ , and hence everywhere. Since it is either cost-reducing or equivalent to ν, it follows that ν cannot be cost-reducing with respect to ν θ , i.e., that ν θ is a cost minimizer with respect to local modifications.

Lemma 4.14. Any c [F0,F1] -locally optimal normal transport ν with rotation number θ = θ ν is equivalent to the shift ν θ .

Proof.

Let v ′′ i = -v ′ i = i for i = 1, 2, . . . . All local modifications ν i = ν v ′ i ,v ′′
i of ν constructed as in Lemma 4.10 cannot be cost-reducing and are therefore equivalent to ν. On the other hand, this sequence stabilizes to the shift ν θ on any bounded subset of V 0 × V 1 as soon as this set is covered by (-i, i) × (-i, i). Therefore ν θ has the same conjugate as all ν i and is equivalent to ν. Lemmas 4.13, 4.14, and 4.5 together imply Theorem 4.6.

5.

Transport optimization for periodic measures. Let now c be a cost function that satisfies the Monge condition (3.1), the integer shift invariance condition (3.2), the growth condition (3.3), and is bounded from below. Suppose that γ θ is a locally optimal transport plan on U 0 × U 1 with marginals µ 0 , µ 1 conjugate to the shift ν θ . Define c [F0,F1] as in (4.2) and let F θ 1 (u) = F 1 (u) -θ as illustrated in fig. 2.1. Definition 5.1. We call the quantity

C [F0,F1] (θ) = 1 0 c [F0,F1] (v ′ , v ′ + θ) dv ′ = 1 0 c F -1 0 (v ′ ), (F θ 1 ) -1 (v ′ ) dv ′ (5.1)
the average cost (per period) of the transport plan γ θ . Observe that it is indifferent whether to integrate here from 0 to 1 or from v to v + 1 for any real v. Examples of average cost functions C [F0,F1] for different marginals µ 0 , µ 1 and different cost functions c are shown in fig. 5.1.

The following technical lemma provides a "bracket" for the global minimum of C [F0,F1] and estimates of its derivatives independent of µ 0 , µ 1 .

Lemma 5.2. The average cost C [F0,F1] is a convex function that satisfies the inequalities

inf x,y c(x, y) ≤ C(θ) ≤ C [F0,F1] (θ) ≤ C(θ) (5.2)
with

C(θ) = inf -1≤u1≤2 θ-1≤u2≤θ+2 c(u 1 , u 2 ), C(θ) = sup -1≤u1≤2 θ-1≤u2≤θ+2 c(u 1 , u 2 ).
(

5.3)

There exist constants Θ < Θ and L, L > 0 such that the global minimum of C is achieved on the interval [Θ, Θ] and

-L ≤ C ′ [F0,F1] (Θ -0) ≤ 0 ≤ C ′ [F0,F1] (Θ + 0) ≤ L, (5.4)
where

C ′ [F0,F1] (•) is the derivative of C [F0,F1]
. These constants are independent on µ 0 , µ 1 and are given explicitly by formulas (5.5), (5.6) and (5.7) below.

The bounds given in the present lemma are rather loose. E.g., for c(x, y) = |x-y| α with α > 1, they are C(θ) = conv min(|θ +3| α , |θ -3| α ), C(θ) = max(|θ +3| α , |θ -3| α ), and -Θ = Θ = 6. For symmetric costs like this one it is often possible to replace [Θ, Θ] by the interval [-1, 1] which may be tighter.

Proof. To prove convexity of

C [F0,F1] it is sufficient to show that C [F0,F1] 1 2 (θ ′ + θ ′′ ) ≤ 1 2 C [F0,F1] (θ ′ ) + C [F0,F1] (θ ′′ ) for all θ ′ , θ ′′ . Let θ ′ < θ ′′ , denote θ = 1 2 (θ ′ + θ ′′ ) and write C [F0,F1] (θ) = 1 0 c [F0,F1] (v, v + θ) dv = θ-θ ′ +1 θ-θ ′ c [F0,F1] (v ′ , v ′ + θ) dv ′ , C [F0,F1] (θ ′ ) = θ-θ ′ +1 θ-θ ′ c [F0,F1] (v ′ , v ′ + θ ′ ) dv ′ , C [F0,F1] (θ ′′ ) = 1 0 c [F0,F1] (v, v + θ ′′ ) dv.
Making the change of variables v ′ = v +θ -θ ′ and taking into account that θ -θ ′ +θ = 2θ -θ ′ = θ ′′ , we get Since v + θ -θ ′ > v and v + θ ′′ > v + θ, the Monge condition for c implies that the integrand here is negative on a set of nonzero measure, yielding the desired inequality for the function C [F0,F1] . Note that convexity of C [F0,F1] implies its continuity because C [F0,F1] is finite everywhere. Bounds (5.2) on C [F0,F1] (θ) follow from (5.1)

2C [F0,F1] (θ) -C [F0,F1] (θ ′ ) -C [F0,F1] (θ ′′ ) = 1 0 c [F0,F1] (v, v + θ) + c [F0,F1] (v + θ -θ ′ , v + θ ′′ ) -c [F0,F1] (v + θ -θ ′ , v + θ) -c [F0,F1] (v, v + θ ′′ ) dv.
with v = 0 because v ′ -1 ≤ F -1 0 (v ′ ) ≤ v ′ + 1, v ′ + θ -1 ≤ (F θ 1 ) -1 (v ′ ) ≤ v ′ + θ + 1
, and 0 ≤ v ′ ≤ 1. Furthermore, the growth condition (3.3) implies that C(θ) ≥ P as soon as |θ| > R(P ) + 3. Indeed, in this case |u 2 -u 1 | ≥ |θ| -3 ≥ R(P ) and right-hand sides of formulas (5.3) are bounded by P from below. Therefore one can set

Θ = inf{θ : C(θ) = min θ ′ C(θ ′ )} > -∞, Θ = sup{θ : C(θ) = min θ ′ C(θ ′ )} < ∞, (5.
(θ) ≥ C(θ) > min θ ′ C(θ ′ ) ≥ min θ ′ C [F0,F1] (θ ′ ), so θ cannot be- long to arg min θ ′ C [F0,F1] (θ ′ ); a similar conclusion holds if θ > Θ.) It follows that C ′ [F0,F1] (Θ -0) ≤ 0 ≤ C ′ [F0,F1] (Θ + 0). By convexity C ′ [F0,F1] (Θ + 0) ≤ (C [F0,F1] (θ) -C [F0,F1] (Θ)
)/(θ -Θ) for all θ ≥ Θ. The right-hand side of the latter inequality can be estimated from above by

L = inf θ≥Θ C(θ) -C(Θ) θ -Θ . (5.6)
The ratio in the right-hand side takes finite values, so L is finite. This establishes the inequality C ′ [F0,F1] (Θ + 0) ≤ L. The rest of (5.4) is given by a symmetrical argument; in particular

L = inf θ≤Θ C(θ) -C(Θ) Θ -θ . (5.7)
Definition 5.3. A locally optimal transport plan γ θ0 is called globally optimal if θ 0 ∈ arg min θ C [F0,F1] (θ).

We can now reduce minimization of (2.1) on the unit circle to minimization of (2.2) on R, which involves the cost function c rather than ĉ:

Theorem 5.4. The canonical projection π : R → T establishes a bijection between globally optimal transport plans on R × R and transport plans on T × T that minimize (2.1).

Proof. A transport plan γ on T × T minimizes (2.1) if it is a projection of a transport plan on R × R that locally minimizes the transport cost defined by the cost function ĉ(x, y) = min k∈Z c(x, y + k) (see introduction; min here is attained because of the integer shift invariance and growth conditions (3.2), (3.3)).

Denote S = {(x, y) : c(x, y) = ĉ(x, y)} and observe that the support of the globally optimal plan γ θ0 lies within S: indeed, if it did not, there would exist a (nonlocal but periodic) modification of γ θ0 bringing some of the mass of each period to S and thus reducing the average cost. Therefore γ θ0 is locally optimal with respect to the cost ĉ(x, y) and its projection to T × T minimizes (2.1).

Conversely, a minimizing transport plan on T × T can be lifted to R × R in such a way that its support lies inside S (translations of arbitrary pieces of support by integer increments along x and y axes are allowed because they leave ĉ(x, y) invariant). Therefore its average cost per period cannot be less than that of a globally optimal transport plan on R × R.

If

θ -θ < ǫ/L, then compute C [F0,F1] (θ), C [F0,F1] (θ), solve the linear equation C [F0,F1] (θ) + C ′ [F0,F1] (θ + 0)(θ -θ) = C [F0,F1] (θ) + C ′ [F0,F1] (θ -0)(θ -θ) (6.5)
for θ, and stop. 6. Otherwise set θ :

= θ if C ′ [F0,F1] (θ + 0) < 0, or θ := θ if C ′ [F0,F1] (θ -0) > 0. 7. Go to step 2.
It follows from inequalities (5.4) of Lemma 5.2 that the minimizing value of θ belongs to the segment [Θ, Θ]. Therefore at all steps

C ′ [F0,F1] (θ + 0) ≤ 0 ≤ C ′ [F0,F1] (θ -0) (6.6)
and the segment [θ, θ] contains the minimum of C.

Step 5 requires some comments. By convexity,

-L ≤ C ′ [F0,F1] (θ ± 0) ≤ L for all Θ ≤ θ ≤ Θ, i.e., |C ′ [F0,F1] (θ ± 0)| ≤ L at all steps. When θ -θ < ǫ/L, this bound ensures that for any θ ′ in [θ, θ] the minimal value of C is within ǫ/L • L = ǫ from C [F0,F1] (θ ′ ).
If there is a single exceptional value of θ in that interval, then it is located precisely at the solution of (6.5) and must be a minimum of C because of (6.6), so the final value of θ is the exact solution; otherwise it is an approximation with guaranteed accuracy.

The final value of θ will certainly be exact when masses of all atoms are rational numbers having the least common denominator M and ǫ < 1/M . Indeed, in this case any interval [θ, θ] of length ǫ can contain at most one exceptional value of θ.

Since at each iteration the interval [θ, θ] is halved, step 5 will be achieved in O(log 2 ((Θ -Θ)/(ǫ/L))) iterations. By Lemma 6.1 each instance of step 3 (and equation (6.5)) takes O(n 0 + n 1 ) operations. Thus we obtain the following result. In the first experiment the value of ǫ was set to 10 -10 , the algorithm was run 10 times for each pair (n 0 , n 1 ) with 1 ≤ n 0 , n 1 ≤ 100, and the computing times were averaged. In the second experiment (n 0 , n 1 ) was fixed at (10, 10) and the average computing time was similarly computed for different values of ǫ. The averaged computing times for the two experiments are plotted in fig. 6.2. Observe the manifest linear dependence of computing time on n 0 + n 1 and log ǫ. The next figure is a concrete, if not entirely serious, illustration of optimal matching in the case of distributions on the "color circle."

Recall that in the HSL (Hue, Saturation, and Lightness) color model, the color space is represented in cylindrical coordinates. The polar angle corresponds to the hue, or the degree to which a color can be described as similar to or different from other colors (as opposed to difference in saturation or lightness between shades of the same color).

We chose two famous paintings, one Russian and one French, whose highly different coloring is characteristic of the two painters, the Expressionist Filipp Maliavin (1869-1940) and the Symbolist Pierre Puvis de Chavannes (1824-1898). An optimal matching of the hue distributions according to the linear cost c(x, y) = |x -y| 2 was used to substitute hues of the first painting with the corresponding hues of the second one while preserving the original values of saturation and brightness. In spite of the drastic change in coloring, the optimality of matching ensures that warm and cold colors retain their quality and the overall change of aspect does not feel arbitrary or artificial.

7. Related algorithmic work. Fast algorithms for the transportation problem on the circle, with the Euclidean distance |x -y| as a cost, have been proposed in a number of works. Karp and Li [START_REF] Karp | Two special cases of the assignment problem[END_REF] consider an unbalanced matching, where the total mass of the two histograms are not equal and elements of the smaller mass have to be optimally matched to a subset of elements of the larger mass. A balanced optimal matching problem has later been considered independently by Werman et al [START_REF] Werman | Bipartite graph matching for points on a line or a circle[END_REF]; clearly, the balanced problem can always be treated as a particular case of the unbalanced one. In both of these works O(n log n) algorithms are obtained for the case where all points have unit mass.

Aggarwal et al [1] present an algorithm improving Karp and Li's results for an unbalanced transportation problem on the circle with general integer weights and the same cost function |x -y|. They also consider a general cost function c(x, y) that satisfies the Monge condition and an additional condition of bitonicity: for each x, the function c(x, y) is nonincreasing in y for y < y 0 (x) and nondecreasing in y for y > y 0 (x). Note that this rules out the circular case. The second algorithm of [1] is designed for bitonic Monge costs and runs in O(n log M ) time for an unbalanced transportation problem with integer weights on the line, where M is the total weight of the matched mass and n is the number of points in the larger histogram.

The algorithm proposed in the present article only applies to the balanced problem for a Monge cost. However it does not involve bitonicity and is therefore applicable on the circle, where it achieves the same O(n log M ) time as the second algorithm of [1] if all weights are integer multiples of 1/M . Although our theory is developed for the case of costs satisfying a strict inequality in the Monge condition, it can be checked that the discrete algorithm works for the case c(x, y) = |x -y|, which can be treated as a limit of |x -y| λ , λ > 1, as λ → 1 [START_REF] Rabin | Transportation distances on the circle[END_REF].

Finally we note that results of [1] were extended in a different direction by Mc-Cann [START_REF] Mccann | Exact solutions to the transportation problem on the line[END_REF], who provides, again in the balanced setting, a generalization of their first algorithm to the case of a general cost of the concave type on the open line. This case is opposite to Monge costs and requires completely different tools. Indeed, for a strictly concave cost such as c(x, y) = |x -y| the notion of locally optimal transport plan on the universal cover does not make sense: concave costs favor long-haul transport over local rearrangements, destroying local finiteness.

. 6 )

 6 Definitions (3.4) and (3.6) mean that F µ , F -1 µ

Fig. 5 . 1 .

 51 Fig. 5.1. Average cost functions C [F 0 ,F 1 ] for some atomic marginals and costs. The cost |x -y| is regarded as lim |x -y| λ , λ > 1, as λ → 1.

Theorem 6 . 2 .j = 1 .

 621 The above binary search algorithm takes O((n 0 + n 1 ) log(1/ǫ)) comparisons and evaluations of c(x, y) to terminate. The final value of θ is within ǫ/L from the global minimum, and C [F0,F1] (θ) ≤ min θ C [F0,F1] (θ) + ǫ. When all masses m the least common denominator M , initializing the algorithm with ǫ = 1/2M leads to an exact solution in O((n 0 + n 1 ) log M ) operations.6.3. Experiments.We tested experimentally the estimates of Theorem 6.2 for time complexity as a function of parameters of the problem. The average computing time of the algorithm for different values of n 0 + n 1 and ǫ is illustrated in fig.6.2. These results have been obtained using the following procedure. For each value of n 0 and n 1 , points {x 1 , x2 , . . . , xn0 } and {ŷ 1 , ŷ2 , . . . , ŷn1 }, which constitute the support of distributions µ 0 and µ 1 , are drawn independently from the uniform distribution on [0, 1] and sorted. The masses m (0) i , i = 1 . . . n 0 and m (1) j , j = 1 . . . n 1 are then drawn from the uniform distribution and normalized such that 1≤i≤n0 m Finally the transport cost is minimized for c(x, y) = |x -y|. The code used to produce this figure is available online at the web site of the OTARIE project http://www.mccme.ru/ ∼ ansobol/otarie/software.html.

  (a) Average computing time, sec., vs n 0 + n 1 for ǫ = 10 -10 . (b) Average computing time, sec., vs log 10 ǫ for n 0 = n 1 = 10.

Fig. 6 . 2 .

 62 Fig. 6.2. Average computing time of the algorithm for different values of (n 0 + n 1 ) and log 10 ǫ.The experiment was performed on a PC with a 3.00 GHz processor.

  (c) Puvis' hues transported to Maliavin's canvas.

Fig. 6 . 3 .

 63 Fig. 6.3. Optimal matching of the hue component of color.

  5) where min is attained because C is continuous. The set arg min θ ′ C [F0,F1] (θ ′ ) lies on the segment [Θ, Θ]. (Indeed, if e.g. θ < Θ, then C [F0,F1]

In what follows the words 'locally finite' defining a transport plan will often be dropped.

* Supported by the French National Research Agency (project ANR-07-BLAN-0235 OTARIE, http://www.mccme.ru/~ansobol/otarie/). The hospitality of UMR 6202 CNRS "Laboratoire Cassiopée" (Observatoire de la Côte d'Azur) is gratefully acknowledged. † LTCI CNRS, Télécom ParisTech. -Dauphine and CEREMADE. § Institute for Information Transmission Problems and UMI 2615 CNRS "Laboratoire J.-V. Poncelet," Moscow; partially supported by the Ministry of National Education of France and by the Russian Fund for Basic Research (project RFBR 07-01-92217-CNRSL-a).

6.

6. Fast global transport optimization. In a typical application, such as the image processing problem described in the introduction, measures µ 0 and µ 1 come in the form of histograms, i.e., discrete distributions supported on subsets X = {x 1 , x2 , . . . , xn0 } and Ŷ = {ŷ 1 , ŷ2 , . . . , ŷn1 } of the unit circle. These two sets may coincide. In what follows we replace X and Ŷ with their lifts to the universal cover X and Y and assume that the points of the latter pair of sets are sorted and numbered in an increasing order:

Denote masses of these points by µ

j ; these are assumed to be arbitrary positive real numbers satisfying 1≤i≤n0 m

(1) j = 1. 6.1. Computation of the average cost and its derivative. Define j(θ) as the index of min{y j : F θ 1 (y j ) > 0} and denote y θ 1 = y j(θ) , y θ 2 = y j(θ)+1 , . . . , y θ n1 = y j(θ)+n1-1 . All the values

belong to the segment (0, 1]. We now sort these values into an increasing sequence, denote its elements by

are uniquely defined and belong to X, Y . It is now easy to write an expression for the function C [F0,F1] :

Observe that, as the parameter θ increases by ∆θ, those v (k) that correspond to values F θ 1 decrease by the same increment. Let F θ 1 (y j0 ) be such a value. As it appears in (6.2) twice, first as v (k) and then as -v (k-1) in the next term of the sum, it will make two contributions to the derivative C ′

Fig. 6.1. Derivation of expressions (6.3), (6.4) for C ′ [F 0 ,F 1 ] (θ ± 0). Thick lines show fragments of complete graphs of F 0 , F θ 1 corresponding to jth terms in (6.3), (6.4), thin dashed line (bottom) marks the common value of F 0 and F θ 1 . Note that X = F -1 0 (F θ 1 (y j )) is equivalently expressed as inf{x :

Proof. Sorting the n 0 + n 1 values (6.1) into an increasing sequence requires n 0 + n 1 -1 comparisons (one starts with comparing F 0 (x 1 ) and F θ 1 (y θ 1 ) to determine v (1) , and after this each of the remaining values is considered once until there remains only one value, which is assigned to v (n0+n1) with no further comparison). At the same time, pointers to x (k) and y (k) should be stored. After this preliminary stage, to find the values for C [F0,F1] and its one-sided derivatives it suffices to evaluate each of the n 0 + n 1 terms in (6.2) and to take into account the corresponding contribution of plus or minus c(x (k) , y (k) ) to the value of C ′ [F0,F1] (θ), paying attention to whether the value of θ is exceptional or not. All this can again be done in O(n 0 + n 1 ) operations.


Ziyong Wang

Chen Zhang

Jason Wu

Bo Li

Anna Chrostowska

Panaghiotis Karamanis

Shih- Yuan Liu

trans-Hydroalkynylation of Internal 1,3-Enynes Enabled by Cooperative Catalysis

A cooperative catalyst system involving Pd(0)/Senphos complex, tris(pentafluorophenyl)borane, copper bromide, and an amine base, is demonstrated to catalyze trans-hydroalkynylation of internal 1,3-enynes. For the first time, a Lewis acid catalyst is shown to promote the reaction involving the emerging outer-sphere oxidative reaction step. The resulting cross-conjugated dieneynes are versatile synthons for organic synthesis and their characterization reveals distinct photophysical properties depending on the positioning of the donor/acceptor substituents along the conjugation path.

Hydroalkynylation, i.e., the net addition of terminal alkyne to an alkyne unit, is an atom-economical and straightforward approach to 1,3-enynes and related conjugated molecules which are widely exploited in bioactive molecules, 1 functional materials 2 and complexity-building transformations. [START_REF] Bao | 2-Activated 1,3-enynes in enantioselective synthesis[END_REF] To this end, significant advances have been achieved in terminal alkyne homo-and cross coupling 4,5 (Scheme 1a). However, the hydroalkynylation of an internal alkyne has been more challenging due to the competing homodimerization, oligomerization and cyclotrimerization of either starting materials as well as regio-and diastereoselectivity issues when unsymmetrical internal alkynes are used. So far, a few stereoselective hydroalkynylation 6,4a of internal alkynes (Scheme 1b) have been reported, with trans-selective examples being particularly rare (Scheme 1c). In 2014, Zhu et al. reported the first example of Pd-catalyzed trans-addition terminal alkyne to N-sulfonyl ynamides. [START_REF] Liu | Palladium-Catalyzed Cross Addition of Terminal Alkynes to Aryl Ynamides: An Unusual trans-Hydroalkynylation Reaction[END_REF] Recently, Lassaletta and coworkers developed a catalyst-counteranion-dependent Aucatalyzed divergent trans-hydroalkynylation of haloalkynes or cis-bromoalkynylation of terminal alkynes. [START_REF] García-Fernández | AuI-Catalyzed Hydroalkynylation of Haloalkynes[END_REF] Another elegant work from Fürstner demonstrated that [Cp*RuCl]4 catalyzed the trans-addition of i Pr3SiCCX (X =H, Cl) across internal alkynes. [START_REF] Barsu | Ruthenium-Catalyzed trans-Hydroalkynylation and trans-Chloroalkynylation of Internal Alkynes[END_REF] Nevertheless, a broadly applicable transhydroalkynylation of internal alkynes with wide scope remains underdeveloped.

In our program to expand the chemical space in ligand design for applications in homogeneous catalysis through BN/CC isosterism, [START_REF] Giustra | The State of the Art in Azaborine Chemistry: New Synthetic Methods and Applications[END_REF] we recently described trans-hydroboration, [START_REF] Xu | Site-Selective and Stereoselective trans-Hydroboration of 1,3-Enynes Catalyzed by 1,4-Azaborine-Based Phosphine-Pd Complex[END_REF] trans-cyanoboration [START_REF] Zhang | Pd-Senphos Catalyzed trans-Selective Cyanoboration of 1,3-Enynes[END_REF] and cis-carboboration of 1,3-enynes [START_REF] Wang | C-Boron Enolates Enable Palladium Catalyzed Carboboration of Internal 1,3-Enynes[END_REF] catalyzed by a 1,4-azaborine biaryl phosphine(Senphos)-Pd complex. Mechanistic studies 14 support an unusual outersphere oxidative addition mechanism [START_REF]Palladium(0)-Lithium Iodide Co-catalyzed Asymmetric Hydroalkylation of Conjugated Enynes with Pronucleophiles Leading to 1,3-Disubstituted Allenes[END_REF] where a Lewis acid (LA) activates a Pd(0)-enyne complex to form a zwitterionic Pd-allyl species (e.g., A in Scheme 1d). 16 Recently, we were able to crystallographically characterize an outer-sphere oxidative addition adduct where the LA is B(C6F5)3. 17 Whereas in our reported borylative transformations the LA activator is employed in stoichiometric amounts as substrates, we hypothesized that a catalytic amount of LA [START_REF] Becica | The roles of Lewis acidic additives in organotransition metal catalysis[END_REF] such as B(C6F5)3 may enable a general and selective addition of nucleophiles and electrophiles across a 1,3-enyne via cooperative catalysis (Scheme 1d). 19 In this communication, we present the successful development of the Pd/Senphos-catalyzed transselective hydroalkynylation of internal 1,3-enynes in the presence of Lewis acid and copper co-catalysts and an amine base (Scheme 1e). The nucleophile in this work is proposed to be the Cu-alkynyl species and the electrophile is a proton. Our catalytic transformation is site-, regio-, and diastereoslective, exhibits broad scope and functional group tolerance, and is an atom-economic and simple method to access the family of cross-conjugated dieneynes. A library of donor/acceptorsubstituted derivatives could be readily prepared and characterized by UV-Vis absorption and emission spectroscopy.

Scheme

Transition Metal-Catalyzed Hydroalkynylation of Alkynes

We selected 1,3-enyne 1 and phenylacetylene as model substrates for our reaction optimization efforts.

We determined that a combination of catalytic amounts of Pd/L1, B(C6F5)3, CuBr, and a stoichiometric amount of TMPH (2,2,6,6tetramethylpiperidine) in toluene at room temperature was crucial to obtain the desired product 2. When any of the catalysts or the TMPH was absent under otherwise identical conditions, the desired hydroalkynylation reaction did not proceed (Table 1, entries 1-6). We then explored the ligand effect on this hydroalkynylation reaction. The use of a sterically more demanding ligand with a larger substituent on the C( 3) position (e.g., R 3 = Et or i Pr) decreases the product yield (entries 7-8). When the R = Cy in Senphos is replaced with R = Ph, a reduction in yield is observed (entry 9 vs. entry 1). A similar detrimental effect of the increasing size of the R 3 substituent is also observed for the arylphosphine (R = Ph) ligand series (entries 9-11). We also determined that the performances of the corresponding carbonaceous ligand CC-L1 and the commercially available structurally related analogue MePhos, [START_REF] Surry | Biaryl Phosphane Ligands in Palladium-Catalyzed Amination[END_REF] are significantly inferior to Senphos L1 (entries 12-13 vs. entry 1), highlighting the distinct electronic structure features of the Senphos ligands enabling this transformation. [START_REF] Xu | Monobenzofused 1,4-Azaborines: Synthesis, Characterization, and Discovery of a Unique Coordination Mode[END_REF] The hydroalkynylation proved to be quite general; a broad range of terminal alkynes serve as suitable substrates, affording the corresponding trans-addition product exclusively (Table 2). Various aryl alkynes (entries 4a-4l) bearing electron-donating or electron-withdrawing groups can readily participate in this reaction to afford product in moderate to good yield. Compatible functional groups include fluoro (entries 4b-d), tertiary amine (entry 4e), methoxy (entry 4f), trifluoromethyl (entry 4g), pinacolborato (entry 4h), ester (entry 4i), nitro (entry 4k) and cyano (entry 4l) groups. In addition, silyl acetylenes (entries 4m and 4n) and terminal aliphatic alkynes (entries 4o-4y) of different carbon chain length containing primary chloride (entry 4r), silyl ether (entry 4t), ketone (entry 4u), ester (entry 4v), phthalimide (entry 4w), nitrile (entry 4x) and estrone moiety (entry 4y) are well accepted. We have obtained the X-ray crystal structure of 4m, thus unambiguously establishing the connectivity and the diastereoselectivity of the reaction product. For ease of purification, a more polar enyne was used as a reactant in some cases (entries 4k, 4l, 4r-4y).

We also evaluated the scope with respect to internal 1,3enynes, and we were pleased to find broad applicability of the optimized reaction conditions and exclusive trans-selectivity of the reaction products (Table 3). For example, the reaction is insensitive to an array of electronically different aryl group on the alkene position R 2 (entries 6a-6g). In addition to aryl groups, the R 2 position also tolerates H (entry 6h) and alkyl substituents (entries 6i-6k), albeit requiring slightly forcing conditions. The reaction is also compatible with various diaryl 1,3-enyne substrates bearing electronically diverse aryl groups at the R 3 position (entries 6l-6r).

The regio-and stereochemistry of the hydroalkynylation of a representative diaryl 1,3-enyne substrate is established by an X-ray structure of the product 6m. In addition to Me and aryl groups at the R 3 position, enynes with longer chain alkyl groups (entry 6s) and those bearing functional groups (entries 6t-6y) are also suitable substrates. We also note that 1,3-enynes derived from

B N R 3 Me PR 2 L1 R 3 = Me, L2 R 3 = Et, L3 R 3 = i Pr, L4 R 3 = Me, L5 R 3 = Et, L6 R 3 = i Pr, R = Cy R = Cy R = Cy R = Ph R = Ph R = Ph Me Me PCy 2

CC-L1

Me

PCy 2

MePhos Senphos

(3) commercial drugs such as ibuprofen, isoxepac, naproxen afford the corresponding products 6w-6y in high yield. The hydroalkynylation products are versatile building blocks in organic synthesis. As shown in Scheme 2, the conjugated diene unit of 6l can participate in a Diels-Alder reaction to produce intermediate 7, which can be oxidized [START_REF] Alsharif | DDQ as a versatile and easily recyclable oxidant: a systematic review[END_REF] to give the arene product 8. [START_REF] Kotha | Design and Synthesis of Benzosultinesulfone as a o-Xylylene Precursor via Cross-enyne Metathesis and Rongalite: Further Expansion to Polycyclics via Regioselective Diels-Alder Reaction[END_REF] In the presence of TBAF, desilylation of 6l proceeds smoothly to afford terminal alkyne 9. Building block 9 cycloisomerizes to the corresponding naphthalene derivative 9a in the presence of PtCl2 catalyst. 24 In the presence of TsN3 compound 9 generates triazole 9b in a click reaction catalyzed by CuTc. [START_REF] Raushel | Efficient Synthesis of 1-Sulfonyl-1,2,3triazoles[END_REF] Furthermore, Sonogashira cross coupling of 9 with aryl and acyl electrophiles furnishes the corresponding products 9c, 9d, 9e, respectively. 26 Scheme 2. Derivatization of Reaction Products.

Table 2. Reaction Scope of Terminal Alkynes. a

a Reaction condition: 1 (0.4 mmol), (COD)Pd(CH2TMS)2 (2 mol%), L1 (2 mol%), for aryl alkynes 3a-3l: 3 (1.0 mmol), B(C6F5)3 (5 mol%), CuBr (5 mol%), TMPH (0.4 mmol), toluene (0.8 mL), 22 C, 12 h; for alkyl alkynes 3m-3y: 3 (0.8 mmol), B(C6F5)3 (10 mol%), CuBr (10 mol%), TMPH (0.8 mmol), toluene (0.6 mL), 22 C, 12 h, unless otherwise noted. Yields of isolated product are reported as an average of two runs. Selectivity is determined by 1 H NMR analysis of purified product. See SI for detailed procedure. b 1.2 mmol 4methoxyphenylacetylene was used. c 24 h reaction time. d 1 (0.4 mmol), (COD)Pd(CH2TMS)2 (5 mol%), L1 (5 mol%), B(C6F5)3 (10 mol%), CuBr (5 mol%), TMPH (0.4 mmol), toluene (0.8 mL), 12 h, 50 C. Table 3. Reaction Scope of Internal 1,3-Enynes. a a Reaction condition: 5 (0.4 mmol), (triisopropylsilyl)acetylene (0.8mmol), (COD)Pd(CH2TMS)2 (5 mol%), L1 (5 mol%), B(C6F5)3 (10 mol%), CuBr (10 mol%), TMPH (0.8 mmol), toluene (0.6 mL), 50 C, 12 h, unless otherwise noted. Yields of isolated product are reported as an average of two runs. Selectivity is determined by 1 H NMR analysis of purified product. See SI for detailed procedure. In preliminary labeling experiments we determined that 51% deuterium is incorporated into the product with deuterated phenylacetylene under standard conditions. With deuterated phenylacetylene in the presence of Et3N (instead of TMPH), we observe 93% deuterium incorporation into the product (Scheme 3). These observations from the deuterium labeling experiments are consistent with the amine serving as a proton shuttle to deliver the proton to the activated enyne.

Scheme 3. Deuterium Labeling Studies.

Scheme 4 illustrates a plausible mechanism for the hydroalkynylation reaction. The active Pd(0)/Senphos complex I is generated from the catalyst precursor (COD)Pd(CH2TMS)2 and the Senphos ligand via reductive elimination, [START_REF] Takahashi | Air-and moisture-stable Xantphos-ligated palladium dialkyl complex as a precatalyst for crosscoupling reactions[END_REF] which is then followed by the coordination of the enyne substrate to furnish complex II. The B(C6F5)3 co-catalyst then activates the Pd(0)/enyne complex II to yield the outersphere oxidative addition adduct III. Complex III subsequently undergoes a protodeboronation [START_REF] Tamke | Frustrated Lewis Pair-Catalyzed Cycloisomerization of 1,5-Enynes via a 5-endo-dig Cyclization/Protodeborylation Sequence[END_REF] Using our methodology, we were able to prepare a library of donor/acceptor substituted cross-conjugated 30 dieneynes 10a-e [START_REF]For the synthesis and reactivity studies of the parental compound, 2ethynyl-l,3-butadiene[END_REF] and characterize their photophysical properties. [START_REF] Tykwinski | For studies of the related cross-conjugated enediynes (Y-enynes)[END_REF] As can be seen from Figure 1, the isomeric compounds 10a-e display a relatively narrower abs range (352 -378 nm) than the em range (435 -481 nm). Compounds 10a, 10b, 10e exhibit a bathochromically shifted em compared to compounds 10c and 10d. Among the investigated compounds, 10a stands out as the molecule featuring the highest quantum yield (Φ = 28%) and as the only structure where its lowest energy transition (abs = 360 nm) does not have the strongest extinction coefficient. Electronic structure and TD-DFT calculations at wB97XD/6-311G(d,p) level [START_REF] Gaussian ; Revisiond | Gaussian Inc[END_REF] (see Supporting Information for details) support the experimental observations.

R 2 R 3 B(C 6 F 5 ) 3 R 1 H CuBr B C P Pd R 3 R 2 B C P Pd B C P Pd • R 3 R 2 B(C 6 F 5 ) 3 B C P Pd • R 3 R 2 H B C P R 3 R 2 H R 1 Pd R 1 BrCu CuBr R 2 R 3 H R 1 H 2 NR 2 HNR 2 B(C 6 F 5 ) 3 HNR 2 H 2 NR 2 R 1 H Palladium Copper I II III IV V
TD-DFT reveals that the lowest-energy absorption (abs) is mostly derived from HOMO-LUMO transition for all compounds. Consistent with experimental data, the transition to the second excited state for 10a is predicted to have a larger oscillator strength than the transition to the first excited state. From the computed electronic structure calculations, the HOMO-LUMO transitions have more charge transfer character for compounds 10a, 10b, 10e whereas for compounds 10c and 10d the first excited state is more consistent with a  to * transition. The distinct electronic structure differences of these frontier orbitals and the stronger charge transfer character may be responsible for the observed bathochromic shift in em for 10a, 10b, 10e in CH2Cl2 vs. 10c and 10d. It is worth noting that 10a is the only structure investigated where the donor and acceptor placement involves a cross-conjugated path. In summary, we have developed a highly selective transhydroalkynylation reaction of internal 1,3-enynes with broad scope and functional group tolerance. In developing the reaction, we demonstrate for the first time that catalytic amounts of Lewis acids such as B(C6F5)3 can be deployed to promote a reaction involving an outer-sphere oxidative reaction step. The resulting cross-conjugated dieneynes serve as versatile building block for a variety of conjugated reaction products. In addition, the facile access to a library of isomeric donor/acceptor substituted dieneynes reveals distinct photophysical properties depending on the positioning of the donor/acceptor substituents along the conjugation path. Overall, the unorthodox reactivity pattern exhibited by this reaction should serve as a general platform to other hydrofunctionalizations of unsaturated compounds, which is a subject of current efforts in our laboratory.

ASSOCIATED CONTENT

Supporting Information

The Supporting Information is available free of charge on the ACS Publications website.

Experimental procedures, compound characterization data, computational and crystallographic information (PDF)

Notes

The authors declare no competing financial interests.

b

  Standard conditions, except: (COD)Pd(CH2TMS)2 (2 mol%), L1 (2 mol%), 22 C. c Standard conditions, except: (COD)Pd(CH2TMS)2 (10 mol%), L1 (10 mol%), B(C6F5)3 (20 mol%), 24 h.

  Scheme 4 illustrates a plausible mechanism for the hydroalkynylation reaction. The active Pd(0)/Senphos complex I is generated from the catalyst precursor (COD)Pd(CH2TMS)2 and the Senphos ligand via reductive elimination,[START_REF] Takahashi | Air-and moisture-stable Xantphos-ligated palladium dialkyl complex as a precatalyst for crosscoupling reactions[END_REF] which is then followed by the coordination of the enyne substrate to furnish complex II. The B(C6F5)3 co-catalyst then activates the Pd(0)/enyne complex II to yield the outersphere oxidative addition adduct III. Complex III subsequently undergoes a protodeboronation[START_REF] Tamke | Frustrated Lewis Pair-Catalyzed Cycloisomerization of 1,5-Enynes via a 5-endo-dig Cyclization/Protodeborylation Sequence[END_REF] by an ammonium cation species (originating from a deprotonative activation of terminal alkyne in the Cu cycle to generate the Cu-alkynyl 29 species) to form IV with concomitant release of the B(C6F5)3 cocatalyst. Transmetalation of the Cu-alkynl species (from Cu cycle) to IV affords intermediate V, which upon reductive elimination delivers the product and regenerates the Pd catalyst I.

Figure

  Figure 1. Normalized Absorption and Emission Spectra and Photophysical Data of 10a-10e in CH2Cl2. a

ACKNOWLEDGMENT

  Research reported in this publication was supported by the National Institute of General Medical Sciences of the National Institutes of Health (NIGMS) under Award Number R01GM136920, the Excellence Initiative of Université de Pau et des Pays de l'Adour I-Site E2S UPPA, and by Boston College start-up funds. We also acknowledge the NIH-S10 (award: 1S10OD026910-01A1) and the NSF-MRI (award: CHE-2117246) for the support of Boston College's NMR facilities. The "Direction du Numérique" of the Université de Pau et des Pays de l'Adour (UPPA), Part of this work was granted access to the HPC resources of [CCRT/CINES/IDRIS] under the allocation 2022-2023 [AD010807031R1] made by GENCI (Grand Equipement National de Calcul Intensif)and Mésocentre de Calcul Intensif Aquitain (MCIA) are acknowledged for the support of computational facilities. Z.W. was supported as a LaMattina Graduate Fellow in Chemical Synthesis.

Table 1 . Control Experiment and Ligand Survey of Internal 1,3-enyne trans-Hydroalkynylation a

 1 Yield and selectivity was determined by 1 H NMR analysis of unpurified reaction mixture with CH2Br2 as the internal standard. Selectivity is defined as the ratio of title product over other isomers.

	Entry	Variation from standard condition	Yield (%) b	Selectivity b
	1	none	88(84) c	>98:2
	2	without [Pd]	<2	N.D.
	3	without L1	<2	N.D.
	4	without B(C6F5)3	<2	N.D.
	5	without CuBr	<2	N.D.
	6	without TMPH	<2	N.D.
	7	L2 instead of L1	74	>98:2
	8	L3 instead of L1	9	>98:2
	9	L4 instead of L1	68	>98:2
	10	L5 instead of L1	61	>98:2
	11	L6 instead of L1	22	>98:2
	12	CC-L1 instead of L1	19	>98:2
	13	MePhos instead of L1	13	>98:2

b c

Yield of isolated product in parenthesis.

  Université de Pau et des Pays de l'Adour, E2S UPPA / CNRS, Institut des Sciences Analytiques et de Physico-Chimie pour l'Environnement et les Matériaux IPREM UMR 5254. Hélioparc, 2 avenue P. Angot, 64053 Pau cedex 09, France Panaghiotis Karamanis -Université de Pau et des Pays de l'Adour, E2S UPPA / CNRS, Institut des Sciences Analytiques et de Physico-Chimie pour l'Environnement et les Matériaux IPREM UMR 5254. Hélioparc, 2 avenue P. Angot, 64053 Pau cedex 09, France

	AUTHOR INFORMATION
	Corresponding Author
	Shih-Yuan Liu -Department of Chemistry, Boston College,
	Chestnut Hill, Massachusetts 02467-3860, United States;
	orcid.org/0000-0003-3148-9147; Email: shihyuan.liu@bc.edu
	Authors
	Ziyong Wang -Department of Chemistry, Boston College,
	Chestnut Hill, Massachusetts 02467-3860, United States;
	orcid.org/0000-0002-6585-8951
	Chen Zhang -Université de Pau et des Pays de l'Adour, E2S
	UPPA / CNRS, Institut des Sciences Analytiques et de Physico-
	Chimie pour l'Environnement et les Matériaux IPREM UMR
	5254. Hélioparc, 2 avenue P. Angot, 64053 Pau cedex 09,
	France
	Jason Wu -Department of Chemistry, Boston College, Chestnut
	Hill, Massachusetts 02467-3860, United States
	Bo Li -Department of Chemistry, Boston College, Chestnut Hill,
	Massachusetts 02467-3860, United States
	Anna Chrostowska -