Etienne Talbot
email: etienne.talbot@univ-rouen.fr

Pierre-Emmanuel E Berche

Denis

D Ledue

Monte Carlo investigation of the magnetic anisotropy in Fe/Dy multilayers

Keywords: Monte Carlo simulation, Heisenberg model, Ferrimagnetic multilayers, Perpendicular magnetic anisotropy PACS: 75.40.Mg, 75.30.Gw, 75.25.+z, 75.50.Gg

Monte Carlo investigation of the magnetic anisotropy in

Fe/Dy multilayers

Introduction

Amorphous multilayers of transition metal -rare earth (TM-RE) compounds display particular interesting properties such as a strong magnetic anisotropy which, in some conditions, may be perpendicular to the layers [1,2]. Up to now, the accurate origin of this anisotropy is not clearly understood and several models have been proposed. These models are based on the anisotropic distribution of the TM-RE pairs [1], dipolar interactions [3], local structural anisotropy [4,5] or single-ion anisotropy [6]. However, no model is able to explain all the experimental features.

This study has been performed by taking into account recent experimental results of the structural and of the magnetic properties of amorphous Fe/Dy multilayers [7,8,9]. In particular, hysteresis loops measured with the magnetic field applied successively in the film plane and perpendicular to it show that the perpendicular magnetic anisotropy is stronger for the multilayers deposited at a higher temperature (600K compared to 320K) [8]. In this temperature range, the structural investigations have shown that the multilayer is alternatively composed of layers with a majority of amorphous Fe and of Fe-Dy alloyed layers (Figure 1) [8]. For T S = 320K, the multilayer is made up of roughly pure amorphous Fe and Dy layers with thin interfaces. Above 700K, the layered structure is destroyed, replaced by mixed grains of α-Fe and Fe 2 -Dy which leads to a planar anisotropy. Altogether, these results show that the perpendicular magnetic anisotropy seems to be correlated in these samples with the existence of Fe-Dy amorphous alloy in the multilayers.

In this work, we have performed Monte Carlo (MC) simulations on amorphous Fe/Dy multilayers in order to investigate the influence of the concentration profile on the magnetic properties. We have considered a model including perpendicular magnetic anisotropy (on average) on a small fraction of Dy sites, the other Dy sites exhibit random anisotropy. Two different concentration profiles have been studied: an abrupt profile that is the multilayer is made up of pure Fe layers and pure Dy layers (close to the experimental profile for T S = 320K) and an experimental profile obtained from tomographic atom probe analysis on a multilayer deposited at 570K (Figure 1). Our goal is to investigate by means of MC simulations the influence of the concentration profile on the global magnetisation, on the global susceptibility, on the sublattice magnetisations and on the magnetisation per plane. Since the strength of the Dy single site anisotropy is not clearly elucidated, we have also investigated the influence of the Dy anisotropy constant magnitude. Moreover, the fraction of Dy sites with perpendicular magnetic anisotropy is not known, so we have also let vary this parameter in our simulations to determine its influence. The model is presented in section 2 and the simulation method is described in section 3. The results are discussed in section 4 and a conclusion is given in section 5.

Model

Magnetic parameters

Our model of Fe/Dy ferrimagnetic multilayers consists of classical Heisenberg spins located at the sites of a face centered cubic lattice with nearest-neighbour exchange interactions. We have chosen this lattice in order to reproduce the experimental observed compact structure [10,11,12]. As the multilayers are amorphous, the exchange interactions are modulated by a gaussian distribution (Figure 2) which models the distribution of the interatomic distances in real samples [13]. On the other hand, the J Fe-Fe and J Fe-Dy exchange interaction values have been extracted from those determined by Heiman et al. [10] on Fe-Dy amorphous thin films by mean-field calculations and adjusted by MC simulations in order to get the pure amorphous Fe Curie temperature (270K [10]). So, they linearly depend on the local concentration X Fe following:

J Fe-Fe (X Fe )/k B = 77 + 449(1 -X Fe ) (in K) (X Fe > 0.4), J Fe-Dy (X Fe )/k B = 8 -198(1 -X Fe ) (in K) (X Fe > 0.4).
Concerning the J Dy-Dy exchange interaction which is known to be much smaller than the others, this latter has been taken constant in function of the concentration. Its value has been estimated by MC simulations in order to provide the pure amorphous Dy Curie temperature (110K [14]): J Dy-Dy /k B = 6.5 (K). These concentration rules are valid only when X Fe > 0.4. J Fe-Dy is strongly negative; it is thus responsible for the ferrimagnetic order (with antiparallel Fe and Dy moments) at low temperature which has been experimentally observed by polarized neutrons reflectivity measurements [9]. We have to mention that the Fe spin value also depends on the local concentration [10]:

S Fe (X Fe ) = 1.1 -1.125(1 -X Fe ) (X Fe > 0.4).
The other numerical values of the model are those of the free ions: g Fe = 2; g Dy = 4/3; S Dy = 2.5. The variation of the exchange energy as a function of the Fe concentration for the three different bonds is shown in Figure 3. We can see that the exchange energy (in absolute value) of the Fe-Fe and Fe-Dy bonds are maximum respectively at X Fe = 0.8 and X Fe = 0.5 and that the contribution of the Fe-Dy pairs will play a significant role in the magnetic ordering.

Concentration profile

The introduction of realistic concentration profiles along the growth direction of the multilayer is essential since atomic diffusion seems to have a major influence on the ) is directly obtained from tomographic atom probe analyses of a multilayer (Fe 3nm/Dy 2nm) built up at 570K [8] (Figure 1). This profile is composed of a rich Fe region (Fe 90 Dy 10 ) and a wide region in which the concentration varies. As it has been previously mentioned, this concentration profile leads to the maximum of perpendicular magnetic anisotropy experimentally observed by Tamion et al. [8,9]. In the follow- ing, this profile will be called TAP profile (for tomographic atom probe profile). The main difference between the two profiles is that the pure Dy region in the multilayer with an abrupt profile is replaced by an alloyed region with variable local concentration in the case of the TAP profile.

Magnetic anisotropy

In order to investigate the perpendicular magnetic anisotropy in the Fe/Dy multilayers, we have focused our attention on a local structural anisotropy model [4] which seems to be well adapted to understand the referenced experimental results.

This model is based on the central hypothesis that during layer-by-layer growth small Fe-Dy crystallized zones at a scale of a few interatomic distances are formed defining on average a preferential axis perpendicular to the film plane. Indeed, they probably exhibit a structure that is reminiscent of those of the defined TM-RE compounds, whose hexagonal symmetry leads to anisotropy axes on average perpendicular to the layers. These locally ordered clusters consist in our model of 13 atoms (one central Fe atom and its 12 nearest neighbours) (Figure 5). Among the 8 neighbours which are not in the plane of the central atom, between 2 and 4 Dy atoms are randomly distributed in order to obtain a cluster concentration close to that of defined compounds. Since RE atoms are characterized by a strong intrinsic anisotropy due to the particular shape of their 4f orbitals, we have thus considered in the Hamiltonian a single site anisotropy term for all Dy atoms. The anisotropy axes on the Dy atoms belonging to the clusters are parallel to the Fe-Dy bonds (Fe here is the central atom of the cluster); so the unit vectors z i of these 4 axes are (±

1 √ 2 ,0, 1 √ 2 ) or (0, ± 1 √ 2 , 1 √
2 ). This leads on average to an anisotropy direction perpendicular to the plane of the layers. For the other Dy sites which do not belong to the clusters we have considered random anisotropy because of the amorphous structure of the multilayers. The energy of the system is thus:

E = - <i,j> J ij (S i .S j ) -D Dy i / ∈clusters (S i .n i ) 2 + i∈clusters (S i .z i ) 2 , (1)
where S i is the spin of the site i, J ij is the exchange interaction between the nearest-neighbour sites i and j, D Dy is the anisotropy constant on the Dy sites, n i is an unit vector indicating the random anisotropy direction on the Dy site i and z i is an unit vector along the local anisotropy axis of a Dy site belonging to the clusters. We also define the magnetic susceptibility:

χ(T ) = m 2 T -m 2 T k B T × N,
where m T is the magnetisation defined as

m T = 1 N N i=1 m x i 2 + N i=1 m y i 2 + N i=1 m z i 2 1/2 T .
As the anisotropy coefficient D Dy in amorphous multilayers is not accurately determined in the literature, we will consider it as a free parameter with the same value for all the Dy sites (in the clusters and in the amorphous matrix). In the same way, the cluster concentration, which is defined as the number of atoms included in the clusters divided by the total number of atoms, is a free parameter in the simulations since it cannot be evaluated experimentally. But in this work, most of the results have been obtained with a quite small cluster concentration of 5%.

Numerical method

The numerical method is the importance sampling MC method in the canonical ensemble [START_REF] Binder | Monte Carlo Simulation in Statistical Physics[END_REF][START_REF] Heermann | Computer Simulation Methods in Theoretical Physics[END_REF][START_REF] Mac Keown | Stochastic Simulation in Physics[END_REF]. The thermodynamic equilibrium at each temperature, corresponding to the minimization of the free energy of the system, is obtained by the Metropolis algorithm [START_REF] Metropolis | [END_REF]. In this algorithm, the spins are examined individually. A site i is randomly chosen and a unit vector u i defined by the random choice with uniform distribution of its z-component u i ∈ [-1, 1] and an its azimuthal angle ϕ i ∈ [0, 2π[ is determined. The new spin S ′ i is then given by:

S x i ′ = S i 1 -u 2 i cos ϕ i , S y i ′ = S i 1 -u 2 i sin ϕ i , S z i ′ = S i u i ,
where S i = S Fe or S Dy . It has to be noted that this spin trial rotation procedure is isotropic. Then the energy variation ∆E associated to this rotation is calculated. The next step is the following:

-if ∆E ≤ 0, the rotation is accepted; -if ∆E > 0, the rotation may be accepted with a probability that is proportional to the Boltzmann factor exp(-∆E/k B T ) in order to take into account the thermal fluctuations. One MC step consists in examining all spins of the system once. At each temperature, 5000 MC steps were performed to reach the thermodynamic equilibrium, and afterwards the physical quantities were measured by averaging over the next 10000 MC steps. As we are interested in disordered systems, it is necessary to perform simulations on several disorder configurations (typically between 20 and 100) in order to estimate the thermodynamic quantities of a macroscopic system. It is important to note that for saving CPU time, we have limited our simulations to a double bilayer. In order to match up to the experimental system (Fe 3nm/Dy 2nm), we have considered 12 Fe and 6 Dy atomic planes in each bilayer because the ratio of the atomic Fe and Dy radii is around 1.5. Each plane contains 800 atoms (which gives a total number of atoms equal to 28800). Since we are not interested in the critical behaviour, finite size effects have no significant influence on the physical properties under consideration.

Numerical results

We present here the magnetic properties of Fe/Dy multilayers with the concentration profiles defined previously (called abrupt and TAP profiles, see Figure 4). As a first step, since we are interested in the influence of the concentration profile on the magnetic ordering, we have considered a multilayer without any anisotropy term. The next step of this study is then focused on the effect of the magnetic anisotropy on the Dy atoms in the clusters and in the amorphous matrix. No anisotropy term has been introduced into the Hamiltonian, this latter is then invariant under any rotation of all the spins. In these samples, the ground states consist of a ferromagnetic Fe layer and a ferromagnetic Dy layer which are antiparallel. These ground states, so called ferrimagnetic, will be considered as a reference in the following. The thermal variation of the magnetic susceptibility for the whole sample and for each sublattice is shown in Figure 6 with the abrupt and TAP concentration profiles. We observe that the abrupt profile is characterized by two separate magnetic ordering at respectively 285K and 150K corresponding to the Fe layer and the Dy layer. The Dy susceptibility curve (Figure 6(b)) exhibits a shoulder at T ∼ 280K related to the polarisation, via the exchange interaction, of the Dy atoms located in the neighbourhood of the interface. For the TAP profile, the two sublattices order at the same temperature (325K). The shoulder of the Dy susceptibility curve below the peak temperature corresponds to the more concentrated Dy planes. The thermal variation of the magnetisation of the whole Fe/Dy multilayer (Fig- ure 7(a)) is strongly influenced by the concentration profile even if the concentration of each specie is the same in the two cases (X Fe = 2/3).

We note in the case of the abrupt profile the presence of a compensation point at T = 175K when the two sublattice magnetisations are antiparallel with the same modulus. Concerning the TAP profile, the magnetic ordering is continuous and smooth because of the mixing of the two species.

Influence of the magnetic anisotropy: sperimagnetic order

Metastable states

The ground states result from the competition between exchange interactions and magnetic anisotropy. The exchange interactions lead to a collinear order in these systems since there is no frustration whereas the magnetic anisotropy favours a angular distribution of the magnetic moments. Then the ground states are expected to be sperimagnetic with a cone angle that should increase with the anisotropy constant D Dy .

As it has been described in the section 2.3, due to the locally ordered clusters, the average orientation of the moments should be perpendicular to the layers. Because of the competition between exchange and anisotropy ener- gies, numerical convergence towards metastable states at low temperature can occur during the simulation, mainly in the case of small cluster concentrations and small values of D Dy . Figure 8 shows the histograms over 100 chemical disorder configurations of the perpendicular magnetisation component of the Dy sublattice for a multilayer with the TAP profile at T = 1K. The concentration of ordered clusters is of 10% and D Dy /k B = 10K. Without any random magnetic anisotropy (RMA) on the Dy atoms of the matrix (Figure 8 In particular, some disorder configurations exhibit a magnetisation which is clearly non perpendicular to the layers. This phenomenon is less pronounced when the anisotropy constant D Dy increases (Figure 8(c)). Anyway, this problem has been overcome, by decreasing the cooling rate during the simulation and also increasing the number of MC steps at each temperature.

Influence of the concentration profile for different D Dy values

Here we have studied the influence of the magnetic anisotropy constant on the Dy sites in a multilayer containing 5% of locally ordered clusters which means that only 3.5% of Dy atoms exhibit uniaxial anisotropy (on average).

• Abrupt profile We firstly presente the results concerning the abrupt concentration profile. Figure 9 shows the thermal variation of the global magnetisation and of the Fe and Dy sublattices for different values of D Dy . We observe a strong influence of D Dy on the low temperature magnetisation (T < T comp = 175K) (Figure 9(a)), especially on the Dy sublattice, its influence on the Fe sublattice being imperceptible (Figure 9(b)). This phenomenon can be understood by the fact that the magnetic ordering process above T comp is essentially due to the Fe magnetic moments which do not display any magnetic anisotropy. Below T comp , both the low-temperature global magnetisation and the Dy sublattice magnetisation decrease when the anisotropy constant increases. This is due to the increase of the non collinearity of the Dy moments, that is a sperimagnetic order already described in the literature in similar systems [19,20]. These magnetic configurations are qualitatively shown in Figure 10 for three values of D Dy . The angular distribution of Dy moment is broader in the core of the Dy layer whereas it is the contrary for the Fe moments which are more misaligned close to the interface because of Fe-Dy coupling.

• TAP profile The influence of the concentration profile on the sperimagnetic order has been investigated, using the previously defined TAP profile. Figure 11 shows the thermal variation of the global magnetisation and of the Fe and Dy sublattices for a concentration of ordered clusters of 5% and different values of D Dy .

As in the case of the abrupt profile, we observe a decrease of the global magnetisation of the multilayer when D Dy increases but this decrease is much less marked than in the previous case. Indeed, we note that the magnetisation of the Dy sublattice is less influenced by D Dy whereas a small effect on the Fe sublattice can be seen when D Dy /k B = 50K. It means that as the layers are more mixed, the influence of the RMA on the Dy atoms is partially transmitted to the Fe sublattice which displays also a low temperature non-collinear magnetic order (Figure 12). We can see that, In order to get a quantitative analysis of the sperimagnetism, we have calculated the magnetisation for each plane and the histogram of the perpendicular component of the magnetic moments for each plane and each sublattice. Since the sperimagnetism is more pronounced for D Dy /k B = 50K, we have chosen this case. The results obtained at T = 1K are shown in Figures 13 and14 for the abrupt and TAP profiles respectively.

For the abrupt profile, we clearly see the evolution of the magnetisation per plane which is not homogeneous inside the multilayer (Figure 13). The decrease of the magnetisation is related to the broadening of the distribution of the perpendicular component of the atomic moments as mentioned previously. More precisely, far from the interface, the distribution is very broad in the Dy atomic planes (Figure 13(a)) due to the large value of the ratio D Dy /J Dy-Dy = 7.7. At the interface, the distribution width is minimum in the Dy layer (Figure 13(b)) whereas it is maximum in the Fe layer (Figure 13(c)) because of the strong Fe-Dy coupling. Finally, the distribution is close to a Dirac peak in the core of the Fe layer since there is no magnetic anisotropy (Figure 13(d)).

In the case of the TAP profile (Figure 14), the variation of the magnetisation per plane for each sublattice is different from the previous case since it depends on the local concentration. We observe that the magnetisation of the Fe sublattice reaches a maximum for the high-concentrated Fe planes (X Fe = 0.9) whereas for the Dy sublattice, the magnetisation per plane exhibits a maximum for a plane with X Fe close to 0.4. This agrees with the distribution of the perpendicular component of the magnetic moments for the Dy sublattice which is significantly sharper in planes with X Fe = 0.4. Contrarily, for the Fe sublattice, the distribution is sharper in the rich Fe planes. In order to explain these features, we have reported in Table 1, the magnitude of each energy term (per atom) defined as:

E Dy RMA = (1 -X Fe -X cluster Dy )D Dy S 2 Dy , E Dy-Dy exch = (1 -X Fe ) 2 J Dy-Dy S 2 Dy , E Fe-Dy exch = X Fe (1 -X Fe )J Fe-Dy S Fe S Dy , E Fe-Fe exch = X 2 Fe J Fe-Fe S 2
Fe , where X cluster Dy is the Dy atomic fraction in the clusters. In the 24 th plane (rich-Fe plane), the coupling between the Dy and Fe sublattices is weak, then the Fe moments are almost collinear whereas the Dy moments are strongly out of line due to the large ratio E Dy RMA /E Dy-Dy exch + E Fe-Dy exch ∼ 5.2. Concerning the 16 th and 13 th atomic planes, the Dy-Fe sublattice coupling is strong and the possible sperimagnetism of the Dy sublattice (due to the very large E Dy RMA value) could be transmitted to the Fe sublattice. The main difference between these two planes is the ratio E Dy RMA /(E Dy-Dy exch + E Fe-Dy exch ) which is about 1.9 and 2.1 in the 13 th and 16 th atomic planes respectively. This explains why the magnetic configuration is more collinear in the 13 th plane than in the 16 th plane. The comparison of these two values with those of the 24 th plane accounts for the more marked sperimagnetism of the Dy sublattice in this latter. These results evidence that the magnetisation profile along the multilayer is not homogeneous, which has already been observed experimentally on these samples from polarized neutrons reflectivity measurements [9]. We would like to emphasize that unlike what is mentionned in several papers [21,22,23], our results indicate that the atomic moment distribution around the mean magnetisation direction in each plane is not uniform but rather gaussian for the two concentration profiles.

Influence of the cluster concentration

The influence of the cluster concentration also has been investigated for D Dy /k B = 50K in the case of the TAP profile. Figure 15 shows the thermal variation of the reduced magnetisation for such a multilayer for different cluster con-centrations. The increase of the cluster concentration from 5% to 20% has no influence on the global magnetisation: the three curves for each concentration are superimposed at each temperature (Figure 15(a)). Looking more precisely at the sublattice magnetisations at low temperature (Figure 15(b)), we note a very small increase of each sublattice magnetisation as the cluster concentration rises. Quite surprisingly, the cluster concentration has no significant effect on the angular distribution of the moments unlike the anisotropy constant to which the magnetic order is very sensitive. This result can be interpreted from the strong influence of the exchange interactions in these multilayers which induce a global behaviour of the magnetic moments as previously seen.

Conclusion

In order to qualitatively describe the macroscopic perpendicular anisotropy, we have used in our model a local structural anisotropy model which induces a perpendicular magnetisation even for very small cluster concentrations. This study has shown using a numerical investigation the influence of the concentration profile on the thermal variation of the magnetisation and on the magnetisation profile for several values of D Dy . For the abrupt profile, a compensation point can be observed and sperimagnetism occurs essentially in the Dy layer. Concerning the TAP profile, the sperimagnetism on the Dy sublattice is not homogeneous and is less pronounced in the rich Dy plane with X Dy ∼ 0.60. In a near future, we plan to perform simulations of hysteresis loops in order to explain experimental data by suggesting rotation mecanisms for the magnetic moments in function of the concentration profile, the anisotropy constant D Dy and the cluster concentration.

Fig. 1 .

 1 Fig. 1. Concentration profile of an (Fe 3nm/ Dy 2nm) multilayer deposited at 570K. The data have been obtained by tomographic atom probe measurements [8].

Fig. 2 .

 2 Fig. 2. Gaussian distribution of the exchange interaction ((a) Fe-Fe exchange interaction; (b) Dy-Dy exchange interaction).

FeFig. 3 .Fig. 4 .

 34 Fig. 3. Exchange energy as a function of the Fe concentration for the three bonds.

Fig. 5 .

 5 Fig. 5. Schematic representation of the ordered Fe-Dy clusters. Fe atoms are in blue, Dy atoms in green. The arrows show the local anisotropy axes on the Dy sites.

4. 1 .

 1 Influence of the concentration profile without anisotropy (D Dy = 0)

Fig. 6 .

 6 Fig.6. Thermal variation of the magnetic susceptibility for a multilayer (Fe 3nm/Dy 2nm) with the 2 studied concentration profiles ((a) for the whole sample, (b) per sublattice). The curves for the TAP profile result from the average over 4 disorder configurations.

Fig. 7 .

 7 Fig. 7. Thermal variation of the reduced magnetisation for a multilayer (Fe 3nm/Dy 2nm) with the 2 concentration profiles ((a) for the whole sample, (b) per sublattice). M ferrimagnetic (a) corresponds to the ground states with antiparallel Fe and Dy layer magnetisations whereas Msat (b) is the saturated magnetisation of each ferromagnetic sublattice.

  (a)), we observe a sharp distribution at |M z Dy /M Dy | ∼ 1 which means that the magnetisation of the Dy sublattice is perpendicular to the layers for all disorder configurations. The introduction of RMA on the Dy atoms of the matrix leads to a broader distribution P (|M z Dy /M Dy |) (Figure 8(b)).

Fig. 8 .Fig. 9 .

 89 Fig. 8. Histograms over 100 disorder configurations at T = 1K of the perpendicular magnetisation component of the Dy sublattice measured with the TAP profile, a concentration of clusters of 10% ((a) without any RMA and D Dy /k B = 10K, (b) with RMA on the matrix Dy atoms and D Dy /k B = 10K, (c) with RMA on the matrix Dy atoms and D Dy /k B = 50K).

Fig. 10 .

 10 Fig. 10. Schematic representation of the magnetic order at low temperature in a multilayer with an abrupt concentration profile for 3 different anisotropy constant values D Dy .

Fig. 11 .Fig. 12 .

 1112 Fig. 11. Thermal variation of the reduced magnetisation for a multilayer (Fe 3nm/Dy 2nm) with a TAP concentration profile containing 5% of clusters for different values of the Dy anisotropy constant D Dy ((a) total magnetisation, (b) Fe and Dy sublattice magnetisations).

Fig. 13 .

 13 Fig. 13. Reduced magnetisation of each atomic plane of a multilayer at low temperature (T = 1K) with the abrupt concentration profile containing 5% of clusters and D Dy /k B = 50K. Distribution of the perpendicular magnetic component in the Dy layer and in the Fe layer ((a) Dy -core plane; (b) Dy -interface plane; (c) Fe -interface plane; (d) Fe -core plane).

Fig. 14 .

 14 Fig. 14. Reduced magnetisation of each atomic plane of a multilayer at low temperature (T = 1K) with the TAP concentration profile containing 5% of clusters and D Dy /k B = 50K. Distribution of the perpendicular magnetic moment component in the Dy sublattice and in the Fe sublattice ((a) plane 16 -X Fe = 0.50; (b) plane 24 -X Fe = 0.90)

Fig. 15 .

 15 Fig. 15. Thermal variation of the reduced magnetisation for a multilayer (Fe 3nm/Dy 2nm) with a TAP concentration profile for D Dy /k B = 50K for different cluster concentrations ((a) total and sublattice magnetisation, (b) magnification of the sublattice magnetisation).

Table 1

 1 Random magnetic energy and exchange interaction energy for 3 typical atomic planes for the TAP profile.

	plane X Fe E Dy RMA (K) E Dy-Dy exch	(K) E Fe-Dy exch	(K) E Fe-Fe exch	(K)
	13 0.65 102.34	4.97	49.25	49.34
	16 0.50 146.88	10.16	61.14	21.77
	24 0.90	29.30	0.41	5.24	86.29

AcknowledgementsThe simulations were performed at the Centre de Ressources Informatiques de Haute Normandie (CRIHAN) under the project No. 2004002. Moreover, the authors are indebted to Dr. Alexandre Tamion for valuable discussions on the experimental results.


Terry Bates

Rhiann Jakubowski

James A Taylor

Evaluation of the Concord Crop Load Response for Current Commercial Production in New York

Keywords: fruit thinning, pruning, Ravaz index, Vitis labruscana, yield management

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

Introduction

Commercial Concord (Vitis labruscana Bailey) grape production in the Lake Erie AVA aims to consistently grow the largest possible crop size that meets an obligate quality standard in any given season. Achieving this goal requires the understanding and management of vine crop load in a cool and relatively short growing region [START_REF] Howell | Sustainable grape productivity and the growth-yield relationship: a review[END_REF]. "'Crop load' is the crop size relative to vine size (estimated as pruning weight or leaf area) and is a measure of the sink:source ratio" (Keller 2010, pg 169). In cool-climate Concord production, the dormant cane pruning mass of one-year-old wood, commonly referred to as pruning weight, is directly related to vine total leaf area (unpublished data derived from Bates 2008, equ. Total Leaf Area (m 2 /m) = 8.31 * Pruning Weight (kg/m) + 5.10, R 2 = 0.90). Pruning weight is often used to represent leaf area in crop load calculations. The most common usage is in the Ravaz Index, which is calculated as the yield to pruning weight ratio (Y:PW) [START_REF] Ravaz | L'effeuellage de la vigne[END_REF], and used to indicate if vines are "subjectively" 3 overcropped, undercropped, or balanced. Crop load has also been measured more directly as the exposed leaf area to fruit weight ratio, rather than Y:PW. This is particularly true in studies that investigate additional canopy management, such as canopy division, with a need to account for changes in exposed, rather than total leaf area [START_REF] Shaulis | Response of concord grapes to light, exposure and geneva double curtain training[END_REF][START_REF] Smart | Sunlight into Wine: A Handbook for Winegrape Canopy Management[END_REF][START_REF] Kliewer | Leaf area/crop weight ratios of grapevines: influence on fruit composition and wine quality[END_REF], Bates 2008). Therefore, Y:PW or leaf area:fruit yield are quantitative measurements of vine crop load; however, what constitutes "vine balance" is qualitative based on region, variety, viticulture production goals, and grape market destination [START_REF] Howell | Sustainable grape productivity and the growth-yield relationship: a review[END_REF][START_REF] Taylor | Considerations on spatial crop load mapping[END_REF].

Crop load management in commercial Concord vineyards is a function of manipulating crop size through practices such as pruning level and fruit thinning or by influencing vine size through water and nutrient availability and uptake. This paper deals primarily with the former by investigating crop size management over a range of vine sizes; however, the importance of the latter in commercial Concord production cannot be understated. Water, nutrient, and rootstock management options have been addressed in other Concord studies and all show how increasing vine size supports larger crop size while maintaining fruit quality due to the importance of vine size in the crop load ratio [START_REF] Pool | Vineyard balance -what is it? Can it be achieved?[END_REF]). Crop size management has been long studied in NY Concord vineyards, and practices have evolved in response to the economic pressures for maximum yield at a minimum acceptable juice soluble solids and the lowest production costs (see Bates and Morris 2009 for summary).

The effect of crop load, measured by Y:PW, on fruit quality and vine size has been characterized with mixed results in Vitis vinifera cultivars. A Carignane crop load study with an 8.0 to 19.6 Y:PW range, under warm climate conditions, showed delayed fruit maturity and 4 reduction in vine size when Y:PW > 10. Below a Y:PW of 10, yield was not a factor in wine quality [START_REF] Bravdo | Effect of crop level on growth, yield, and wine quality of a high yielding Carignane vineyard[END_REF]. Crop load studies with Cabernet sauvignon in Oakville, CA (Y:PW range of 3 to 14) also supported Y:PW between 5-10 as being indicative of vines with balanced vegetative and reproductive growth [START_REF] Kliewer | Leaf area/crop weight ratios of grapevines: influence on fruit composition and wine quality[END_REF]. Crop load values < 5 were considered undercropped, did not respond to additional crop reduction, and were at risk of excessive canopy growth and cluster zone shading. Crop load values > 10 were considered overcropped with delayed fruit maturation and retarded vegetative growth.

Several cool-climate crop load studies have questioned the appropriateness of these crop load indicators in regions where more leaf area per unit fruit may be needed to achieve optimum fruit maturity and vegetative growth [START_REF] Howell | Sustainable grape productivity and the growth-yield relationship: a review[END_REF]; however, these studies have provided little evidence to deviate from the general crop load response. Fruit thinning studies with Pinot noir in Oregon and Merlot in Hawkes Bay, New Zealand (Y:PW ranges 0.5 to 3.7 and 3.4 to 9.0, respectively) indicated excessive vegetative growth and canopy density in the fruiting zone as having a larger impact on fruit composition than crop load and questioned the need for fruit thinning at these values [START_REF] King | Timing of Crop Removal Has Limited Effect on Merlot Grape and Wine Composition[END_REF][START_REF] Reeve | Vineyard floor management and cluster thinning inconsistently affect 'Pinot noir' crop load, berry composition, and wine quality[END_REF]. Riesling crop load studies in NY and WA (Y:PW ranges 2.9 to 9.9 and 8.2 to 11.4, respectively) both showed slight delays in fruit ripeness at the highest crop load levels but questioned if the minor improvement in fruit quality from fruit thinning justified the economic cost in labor and crop loss [START_REF] Keller | Cluster Thinning Effects on Deficit-Irrigated Winegrapes[END_REF][START_REF] Preszler | Cluster thinning reduces the economic sustainability of Riesling production[END_REF]. Furthermore, differences in fruit composition and wine quality in these studies were attributed more to seasonal differences than crop load at Y:PW between 5 and 10. Higher crop loads than those recommended (Y:PW > 12) have been reported for some cool-climate V. vinifera and French-American hybrids without detrimental effects on fruit quality; however, these required canopy division to optimize leaf area exposure [START_REF] Reynolds | Influence of grapevine training systems on vine growth and fruit composition: a review[END_REF]. This is similar for Concord on Geneva Double Curtain training but only on excessively large vines, which is not the norm for commercial Concord vineyards [START_REF] Shaulis | Response of concord grapes to light, exposure and geneva double curtain training[END_REF] In New York Concord production studies, which evaluated different balanced pruning formulas or fruit thinning to manipulate crop size in response to vine size, Y:PW ranged from 5 to 20 [START_REF] Shaulis | The interaction of resistant rootstock to the nitrogen, weed control, pruning, and thinning on the productivity of Concord grapevines[END_REF][START_REF] Shaulis | Responses of the grapevines and the grapes to the spacing of and within canopies[END_REF][START_REF] Poni | Interactions of crop level and late season water stress on growth and physiology of field-grown concord grapevines[END_REF]. Crop loads > 10 had lower juice soluble solids than the maximum, eventually leading to the pruning severity recommendation that maintained crop loads between 7 and 10. These crop loads in NY Concord ensured acceptable fruit maturity in any given season and tended to increase vine size over time. Under current production costs and market value for juice processing Concord, the economic breakeven point for producers is roughly 13.5 t/ha (6 tons/acre) and conservative balanced pruning only reaches those yields at large vine sizes, which is why historical recommendations for managing 'Concord' vine size have targeted PW > 1.2 kg/vine (~0.5 kg/m of trellis, [START_REF] Shaulis | Looking ahead with grape growers[END_REF]). Commercial Concord producers, under current market conditions, consider Y:PW < 10 as uneconomically sustainable because of low yields; however, consistently pushing high crop loads results in low, and possibly unmarketable, fruit quality while reducing vine size and return crop potential.

In a New York Concord study, which measured the effect of pruning severity and yield on the harvest date for 16 Brix fruit, the industry standard, the harvest date was delayed at crop loads Subjectively, vine balance for Lake Erie AVA Concord grown for the juice grape market is defined as growing the largest possible crop reaching 16 Brix by commercial harvest (30 to 40 days after veraison) and having no net change in vine size. Grower experience with crop thinning also indicates a seasonal influence on the thinning level decision. The purpose of this study was to (a) identify the Y:PW in Lake Erie Concord that matches a commercial definition of vine balance, (b) determine if the method of crop size management (pruning or fruit thinning) effects the crop load response on Brix and the seasonal change in PW, and (c) understand how seasonal conditions change the target Y:PW for vine balance. Results will be used to improve seasonal pruning and fruit thinning recommendations and give validation to spatial crop load vineyard maps.

Materials and Methods

This research investigated the response of Concord grapevines to crop load in two different but related studies. The first study compared the crop load response of Concord vines pruned to various levels against vines pruned to a relatively high node number and then fruit thinned mid-7 season over a four-year period (2001)(2002)(2003)(2004). The second study investigated the crop load response of mid-season fruit thinning over an 11-year period (1999-2009). These will be referred to as the "Pruning vs. Thinning" and the "Thinning Response" studies, respectively.

Vineyard Description. These studies were conducted in two neighboring 0.41 ha Concord vineyards planted in 1956 and separated by less than 70 m at Cornell's Taschenberg Vineyard Laboratory in Fredonia, NY (42°27'00"N, 79°18'42"W). The experimental block elevation was 231 m above sea level with a 1 to 2% slope and a southeast aspect. The soil was well drained Chenango gravel-loam (USDA 1994) with a surface soil pH of 5.5 and 2% organic matter.

The Lake Erie grape production region is characterized by cool and humid conditions. The 92-year (1926 to 2018) average GDD accumulation (base 10 o C) from 1 April to 31 October recorded at the experiment site was 1532 GDD (Table 1). Average precipitation accumulation for the same time period was 630 mm, with rainfall being evenly distributed throughout the growing season. Supplemental irrigation was not used in these vineyards because it is rarely used for commercial Concord production in NY. The 39-year average bud break was 3 rd May, and the 53year bloom and veraison phenology dates were 13 th Jun and 22 nd Aug, respectively. The seasonal GDD accumulations (1 April to 31 October, base 10 o C) for the course of the 1999-2009 experiment ranged from 1343-1811. The warm 2005 season was the only year to be more than two standard deviations from the long-term GDD mean. The seasonal precipitation accumulation (1 April to 31 October) over the same period ranged from 390-950 mm. The dry 2007 season was the only year to be more than two standard deviations from the 92-year mean. Bloom and veraison were relatively consistent in 10 of the 11 years and were within six days of the long-term average. The 8 exception to this was in 2003, when bloom and veraison were 11 and 10 days later than the longterm mean, respectively.

The own-rooted Concord vines in both blocks were 43 years old at the beginning of the experiment and planted at a row-by-vine spacing of 2.7 x 2.4 m, with rows oriented in an east-west direction. Vines were cordon-trained to a trellis wire at 1.8 m. Floor, nutrient, pest and disease management were done according to commercial standards for western NY Concord vineyards [START_REF] Jordan | Cultural practices for commercial vineyards[END_REF]. No-till weed management was used by maintaining a 1.2 m-wide weed-free zone under the vines with pre-and post-emergence herbicides and treating row centers with one glyphosate application at bloom. Ammonium nitrate fertilizer was surface broadcasted at a rate of 56 kg per ha of actual N across the block in a single application near bud break. Fungicide and insecticide materials and application rates were done according to the New York and Pennsylvania Pest Management Guidelines for Grapes [START_REF] Weigle | New York and Pennsylvania Pest Management Guidelines for Grapes[END_REF]) and varied annually depending on seasonal weather conditions.

Pruning vs. thinning study. For pruning treatments, two balanced pruning levels, 33 + 11 and 66 + 11, and one 100 node pruning level, were used to achieve a range of crop levels over a range of vine sizes. For balanced pruning, vines were rough pruned to approximately 100 nodes, and the weight of one-year-old canes was recorded using a spring scale. The number of retained nodes was then adjusted to retain 33 or 66 nodes for the first 500g of pruning weight and 11 additional retained nodes for each additional 500g of pruning weight (imperial units = 30 + 10 or 60 + 10 per pound pruning weight). For the 100 node pruning treatment, vines were pruned to 100 retained nodes regardless of vine pruning weight. There were 36 count vines per pruning treatment grouped in six-vine plots in six randomized complete blocks; however, for this crop load response 9 investigation, the individual vine measurements were used and the vine data regrouped for analysis based on crop load values rather than spatial location. The same vines received the same treatments for each of the four years.

For thinning treatments, 100 vines each year were pruned to have 120 retained nodes per vine, regardless of vine pruning weight. This was chosen to achieve a relatively high crop potential across a range of vine size classes. At 30 days after bloom, 25 individual vines were randomly selected from the original 100 vine population each year, and the fruit was removed and weighed.

The mean crop size was used to calculate the weight of fruit needed to be removed to achieve 0, 25, and 50% crop reduction on the remaining 75 vines. On the vines to receive 25 or 50% crop removal, clusters were randomly removed at 30 days after bloom from individual vines and weighed until the target weight was achieved [START_REF] Pool | Predicting and controlling crop on machine or minimal pruned grapevines[END_REF][START_REF] Fendinger | Effect of mechanical thinning minimallypruned "Concord" grapevines on fruit composition[END_REF]. Thinning Data analysis. The four-year relationship between pruning weight and yield at each pruning severity was done by binning observations for each year and pruning severity into groups at an interval of 0.1 kg/vine and then calculating the mean for each year. Linear and exponential regression models were fitted in JMP (JMP Pro v13.1, SAS Institute Inc, 2016). For the four-year Pruning vs Thinning crop load analysis, multiple linear regression was performed using standard least squares to investigate the effect of treatment (pruning or thinning), year, and Y:PW on juice soluble solids and the seasonal change in pruning weight. In both cases, year had an effect (p value <0.0001); therefore, within year effects of treatment and Y:PW were further analyzed. To investigate the crop load response in each year, observations were binned into groups at intervals of Y:PW = 2, and linear regression was performed against harvest juice soluble solids and the seasonal change in pruning weight. The general crop load response was generated by calculating 11 the mean for each bin (n = 4 years) and then performing linear regression analysis on the combined pruning and thinning data. Combining the data was done to understand the crop load response, regardless of whether it is achieved via altering the source (leaf) or sink (fruit) components.

Different trials were plotted with different markers although, as hypothesized, the thinning and pruning treatments showed similar crop load response. For the crop load -juice soluble solids relationship, a segment linear regression (broken-stick) model was fitted. Observations from the 11-year thinning were similarly binned in intervals of Y:PW = 2 and the same linear and brokenstick regression fitting was performed for the general crop load responses (n = 11 years). The linear response equations for each year were used to investigate the effect of seasonal GDD and precipitation on the annual crop load response. Multivariate regression using standard least squares in JMP was used to investigate the effect of GDD and precipitation on the Y:PW to achieve 16 Brix or no seasonal change in pruning weight. To investigate the effect of year-one Y:PW on yeartwo yield, the 11-year thinning observations, also binned in intervals of Y:PW = 2, were further grouped into small (0-0.6 kg/vine), medium (0.6-1.2 kg/vine), and large (1.2-1.8 kg/vine) vine size classes (roughly representing one-, two-, and three-pound vines, a common industry reference on standard commercial spacing) and subject to linear regression.

Results

Pruning weight in the own-rooted Concord vines from 2001-2004 ranged from 0.2 to 1.8 kg/vine across all pruning levels (Figure 1). As expected, crop size increased as vine size increased, but the response was not linear and differed depending on pruning severity. Severe balanced pruning (33 + 11) had the lowest yield across all vine sizes because there were fewer retained 12 nodes at any given vine size compared to the other pruning treatments. Balanced pruning retained more fruiting buds as vine size increased; therefore, shoot number, cluster number, and total fruit yield increased as vine size increased (data not shown). This was observed with 33 + 11 pruning; however, crop size did not fully plateau at the largest vine size, indicating that pruning had limited yield potential in this treatment. Balanced pruning at 66 + 11 had a similar vine size to yield response as the 33 + 11 treatment but with higher yields because of the less severe pruning formula.

The 66 + 11 treatment also did not fully plateau but approached the same yield as the 100 node pruning at the highest vines size.

Less severe pruning at either 100 or 120 nodes/vine had higher yields than either of the balanced pruning treatments across the range of vine size. The shape of the vine size to yield response was slightly different for the fixed node treatments compared to the balanced treatments.

In the fixed node treatments, there was a stronger positive response of yield to vine size at low pruning weight up to approximately 0.6 kg PW/vine. For both fixed node treatments, the response plateaued at approximately 1 kg PW/vine at a yield value of 13.6 kg fruit/vine and 17.6 kg fruit/vine for the 100 and 120 node pruning treatments, respectively.

Pruning vs. thinning study. The range of pruning severity and vine size created a population of Concord vines with crop loads (Y:PW) from 2 to 50 (Figure 2). Severe balanced pruning (33 + 11) tended to have the lowest crop load values, while 120 node pruning tended to have the highest crop load values. In addition to the pruning treatments described, 120 node vines were fruit thinned at multiple levels to adjust the crop load down and increase the crop load range.

For comparison, the 33 + 11, 66 + 11, and 100 node pruned treatments were grouped to give a 13 range of crop loads achieved through pruning severity. The 120 node pruned vines with and without fruit thinning were grouped to give a range of crop loads achieved through fruit thinning.

There was a negative linear relationship between crop load and juice soluble solids for both pruning and thinning treatments in each year, except for thinned vines in 2001 (Figure 2). The mean gradient, excluding 2001 thinned vines, was -0.076, indicating that an increase in crop load by 13 decreased juice soluble solids by 1 Brix at harvest. There was an effect of Y:PW on juice soluble solids in every year, an effect of treatment in 3 of 4 years, and an interaction effect in 2001 and 2003 (Table 2). In 2001, the thinned vines had lower average yield and higher average pruning weight than the other three years, which compressed the crop load values to the lower end of the scale. In 2003, the response of juice soluble solids to crop load was greater for pruned vines than thinned vines. For the general response over four years, there was an effect of Y:PW on juice soluble solids but not a treatment effect nor an interaction effect on juice soluble solids (Table 2).

Similarly, there was a negative linear relationship between crop load and the seasonal change in vine pruning weight (Figure 3). The mean gradient for all years and treatments was -0.014, indicating an increase in crop load by 7.1 decreased vine pruning weight by 0.1 kg. There was an effect of crop load on the change in pruning weight for each year and smaller effect of treatment in 2002 and 2004, where thinned vines had a greater increase in pruning weight (Table 2). There was not a treatment effect nor an interaction effect on delta pruning weight when all years were combined. The crop load response was very similar whether the crop loads were achieved through pruning or through thinning treatments (Figure 4).

Thinning response study. From 1999-2009, crop load was calculated on vines pruned to 120 nodes and fruit thinned to different crop levels. In each season, there was a negative linear relationship between crop load and juice soluble solids and between crop load and the seasonal change in pruning weight (Figure 5). The mean gradient for the crop load to juice soluble solids response was -0.07 and the mean gradient for the crop load to pruning weight change was -0.016, both similar to the four-year pruning vs thinning study.

Using the average linear regression equations from the crop load vs. juice soluble solids and vs. pruning weight change (Figure 5), the Y:PW that achieved 16 Brix or no PW change in each season was calculated. The effect of seasonal growing degree days and precipitation on the 15 predicted Y:PW to achieve 16 Brix or no seasonal change in PW was investigated by multivariate regression. There was an effect of GDD on the Y:PW to achieve 16 Brix, but no effect of precipitation or an interaction effect (Table 3). In contrast, the main effect in the multivariate model to predict the seasonal change in PW was precipitation, and GDD had no direct predictive power for the seasonal change in PW (Table 3). The univariate relationships between Y:PW for 16 Brix, Y:PW for no PW change, GDD and precipitation are plotted in Figure 6. In general, in warmer seasons, 16 Brix could be achieved at a higher crop load (Figure 6A). The two warmest seasons, 2005 and 2007 achieved 16 Brix at Y:PW > 40. In contrast, for the seasonal change in pruning weight, wetter seasons were able to attain no seasonal change in PW at higher crop loads (Figure 6D).

To investigate the effect of crop adjustment on return crop potential, vines in the thinning study were pruned to 120 nodes, fruit thinned in year one, again pruned to 120 nodes, and received no other crop adjustment in year two. Vines from the thinning study were further grouped into small (0-0.6 kg/vine), medium (0.6-1.2 kg/vine), and large (1.2-1.8 kg/vine) vine size classes to investigate the effect of year-1 crop load on year-2 crop yield. Within a vine size class, there was no effect of year-1 crop load on return crop yield over a wide range of crop load values (Figure 7).

Year-1 crop load influenced the change in vine pruning weight in year-1 (Figures 3, 4 and5), but it was the absolute vine pruning weight at the beginning of year 2 which influenced yield potential in year 2.

Discussion

The main objectives of this study were to characterize the crop load response in NY Concord, determine if the crop load response was influenced by pruning or fruit thinning, and understand seasonal impacts on crop load targets. The results from this study support the general Y:PW reported in accepted crop load theory (Bravdo et al. 1985, Kliewer and[START_REF] Kliewer | Leaf area/crop weight ratios of grapevines: influence on fruit composition and wine quality[END_REF], but add details specific for Concord production for the juice grape industry in NY and PA.

Maximum juice soluble solids were achieved at Y:PW of approximately 7.2. Crop loads below 7.2 did not result in an increase in juice soluble solids but did continue to increase vine size. This is consistent with other fruit thinning studies that showed little response of fruit maturation to crop load at Y:PW values below 10 and attribute fruit quality differences in this range to regional climatic conditions or seasonal variation [START_REF] Keller | Cluster Thinning Effects on Deficit-Irrigated Winegrapes[END_REF][START_REF] Preszler | Cluster thinning reduces the economic sustainability of Riesling production[END_REF][START_REF] Reeve | Vineyard floor management and cluster thinning inconsistently affect 'Pinot noir' crop load, berry composition, and wine quality[END_REF]). The crop load -Brix inflection point of 7.2 in this cool-climate study was slightly lower than the 10-12 reported in warmer climate studies [START_REF] Kliewer | Leaf area/crop weight ratios of grapevines: influence on fruit composition and wine quality[END_REF] to 13.0) with juice soluble solids between 16.6 and 17.0 and a seasonal increase in vine size. This pruning management is appropriate for conservative producers who want to ensure a ripe crop in any given season and do not want to adopt fruit thinning as a crop control strategy. Since retained nodes and yield increased with increasing vine size under this treatment, economic revenue is maximized at high vine size > 1.2 kg/vine (0.5 kg/m of canopy).

Fixed node pruning at 100 nodes/vine without additional crop adjustment had a mean crop load of 17 (quartile range = 12.0 to 20.0) and tended to be on the high end with small vines and the low end with large vines. In commercial Concord operations in NY, the term "grower pruning" typically refers to pruning standard spaced vines (2.7 m row x 2.4 m vine) to between 80 nodes on 18 small vines and 100 nodes on large vines. This study indicated that 100-node pruning management reasonably satisfies the crop load goal for achieving 16 Brix fruit while maintaining vine size in NY Concord. The disadvantage of high node number "pruning only" crop management is that the crop potential is set during the dormant season and does not allow for in-season crop load adjustments. In the event of frost or poor fruit set, limiting bud number limits potential yield by limiting the number of secondary shoots, thereby undercropping the vines and limiting revenue.

In contrast, high fruit set combined with a cool growing season, as experienced in 2003, results in overcropped vines with low juice soluble solids and a reduction in vine size. In this case, midseason fruit thinning is an option to bring vines into balance.

Managing vines to a relatively high level of crop stress together with the physical stress of fruit thinning makes high node number pruning a questionable treatment. In two of the four pruning vs thinning comparison years, and in eight of the eleven thinning study years, vines pruned to 120 nodes without additional crop adjustment had Y:PW values > 20 with reductions in juice soluble solids and vine size. In Lake Erie Concord production, there is a trend to retaining relatively high node numbers as a result of mechanized pruning and frost risk mitigation. Mechanized pruning reduces production costs but is less precise in managing retained node quantity and quality; therefore, producers err on the side of retaining too many buds. In this study, fruit thinning one month after bloom reduced Y:PW, increasing harvest juice soluble solids and the seasonal change in pruning weight. The Concord crop load response was the same whether a given Y:PW was achieved through pruning severity or by retaining additional fruiting nodes followed by midseason fruit thinning.

19

In practice, Concord crop estimation and adjustment are done by mechanical fruit thinning with a harvester at approximately one month after bloom. The decision to fruit thin is primarily dictated by an assessment of vine crop load and vineyard goals for balance, as described in this study. Deciding how much fruit to thin, or retain, can be additionally adjusted by seasonal climatic conditions. The thinning data in this study indicated that warmer than average seasons achieved 16 Brix at higher Y:PW than cooler than average seasons. Assuming GDD accumulation onemonth after bloom is reflective of total season GDD, high GDD accumulation at fruit thinning means the vineyard can be managed to higher crop load levels and may require less thinning or no thinning at all. Low GDD accumulation at fruit thinning would warrant more fruit thinning to maintain higher sugar accumulation rates from veraison to harvest.

There was also an effect of seasonal precipitation on the change in vine pruning weight.

Holding everything else equal, it is common to observe higher vine size in wet years and lower vine size in dry years. It raises an interesting question concerning the error associated with using pruning weight to estimate vine capacity. If a vine with a Y:PW of 17.5 has a vine size of 1.0 kg in an average precipitation year, 1.2 kg in a wet year, and 0.8 kg in a dry year, is the vine capacity different, possibly as a result of decreased photosynthetic capacity in dry years, or an artifact of growth habit, such as longer internode lengths in wet years?

On a broader scale, vine size did influence fruit yield potential at all pruning levels, as expected (Figure 2), but the return crop analysis in this study also showed that it did not matter how the vine size was achieved (Figure 7). For example, overcropping large vines in year 1, undercropping small vines in year 1, and balancing medium size vines in year 1, all resulted in medium size vines having the same yield potential in year 2 at a given pruning level. This illustrates 20 the continual management interaction between crop load and the vine size-yield response, two concepts often confused. Crop load, managed as the Y:PW, influences if a vine will gain or lose pruning weight during a season; however, it is the absolute vine size (or vine capacity) during the dormant season which determines the next season's yield potential.

Based on the crop load relationships and seasonal responses in this study, crop load ranges and recommendations were developed for NY Concord vine balance (Table 4). Climate change predictions for the northeast US point toward warmer and longer growing seasons, longer spring frost windows, increased spring precipitation, and more frequent water shortages in the growing season [START_REF] Wolfe | Unique challenges and opportunities for northeastern US crop production in a changing climate[END_REF]). If managed properly, both warmer and wetter conditions indicate an opportunity to produce Concord grapes in the future at higher crop loads while maintaining high juice soluble solids and vine size. The potential increase in spring frost events is a risk which can be partially mitigated by retaining high node numbers, either manually or mechanically, and then adopting mid-season fruit thinning to balance vines if needed.

Conclusions

The response of Concord in New York to crop load, measured as Y:PW, was similar to previously reported studies with respect to harvest juice soluble solids and the seasonal change in vine pruning weight. There was a negative relationship between crop load and juice soluble solids, which reached a maximum at a Y:PW of 7.2, and with no further increase below 7. provides this information to the growers for crop estimation and management through the Lake Erie Regional Grape Program. The Concord crop load model developed here will integrate with this information to assist growers with mid-season crop estimation and fruit thinning management.

The crop load response information from this study can also be used to provide viticulture context to spatial vineyard crop load mapping where NDVI and Yield monitor spatial sensor data can be used to calculate spatial vineyard crop load maps.

Literature Cited

Bates T, Morris J. ). Pruning and thinning data were combined for the mean response. A broken stick model was fitted to illustrate the inflection point on the juice soluble solids response (see Equation 1). For the pruning weight response, a linear model was used (see Equation 2). Dashed lines indicate 95 % mean confidence limits.

>

  11 but only by one week at crop loads as high as 20 (Bates 2008). Arguably, this is an acceptable harvest delay for higher crop returns. In a Concord mechanization study, which investigated mechanical shoot and fruit thinning practices (Y:PW range 4 to 27), Y:PW values as high as 22 produced fruit with 16 Brix but reduced vine size. However, crop loads between 11 and 14 have 6 produced acceptable juice soluble solids and maintained vine size (Bates 2017). Recent spatial crop load mapping and variable-rate fruit thinning studies in commercial NY Concord observed Y:PW from 15 to 45 in unthinned blocks. Mechanical fruit thinning reduced mean Y:PW in separate management classifications from 30.1 to 19.8 and produced fruit above 16 Brix for the whole vineyard block (Bates et al. 2018, Taylor et al. 2019). The high crop load values, considered as substantial overcropping in other studies, and the wide range in crop load values under commercial Concord production illustrates the impact of market pressures on vineyard management.

  treatments were assigned randomly to the vine population.Different from the pruning treatments, the thinning treatments were applied on a new vine population each year. For example, vines pruned to 120 nodes and fruit thinned in 2001 were again pruned to 120 nodes but not thinned in 2002 to investigate the effect of altering crop load by fruit thinning in year one on vine growth and the return crop size in year two. A new set of vines, which were untreated in 2001, were pruned to 120 nodes and subjected to the thinning treatments in 2002. This process continued on a new vine population for each of the 11 experiment years. Thinning response study. The 120 nodes/vine pruning with mid-season manual fruit (random cluster) thinning procedure described above was applied from 1999 to 2009. The number of vines used each year varied from season to season based on labor availability. From 1999 to 2004 and then again in 2006 and 2010, 100 vines (25 per crop level treatment) were used. In 2005, 10 the vine number dropped to 24 (six vines per crop level treatment). From 2007 to 2009, 56 vines (14 vines per crop level treatment) were used in each year. Vine measurements. Dormant vine pruning weight was measured on individual vines each year. Harvest date was determined when the juice soluble solids in the lower yielding treatments were between 16 and 17 Brix (Table 1). On average, this occurred at approximately 40 days after veraison and coincided with the first or second week of the commercial Concord harvest in the region. A 100-berry pre-harvest sample was randomly collected from two vertical planes transecting the east and west side of each vine. Juice soluble solids were measured with a handheld refractometer (Leica model no.10423, Leica, Inc., Buffalo, NY), and the clusters from each vine were manually harvested and weighed. Crop load was calculated as vine fruit weight divided by vine pruning weight at the end of a season (Y:PW).

  The four-year pruning and thinning data were combined to characterize a general Concord crop load response. A segmented regression (broken-stick) model was fitted to the Y:PW by harvest juice soluble solids relationship to capture the plateau at the lower Y:PW values (Figure 14 4A) as described by Equation 1. The general response indicated that maximum possible juice soluble solids occurred at or below a Y:PW of 7.2. The crop load to reach the industry standard of 16 Brix, however, occurred at a Y:PW of 20, on average. 𝐼𝐼𝐼𝐼 Y: PW < 7.2 𝑡𝑡ℎ𝑒𝑒𝑒𝑒 Brix = 17.2 𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑒 Brix = 17.8 -0.085 * (Y: PW) � Equation 1 R 2 = 0.92, F ratio = 384.8, Prob > F; <0.0001 The Y:PW by seasonal change in vine pruning weight did not plateau but exhibited a linear response across the full range of crop load values (Figure 4B) as described by Equation 2. On average, a Y:PW of 17.5 resulted in no net seasonal change in vine pruning weight. ∆PW(𝑘𝑘𝑘𝑘) = 0.21 -0.012 * (Y: PW) Equation 2 R 2 = 0.68, F ratio = 70.97, Prob > F; <0.0001

  , supporting the theory that cool-climate grapevines need relatively more leaf area per unit fruit to reach maximum maturity. Given the wide range of Y:PW values, this slight difference does not practically change the accepted industry standard crop load range of 5 -10 for high fruit maturity.This study also supports the theory that vines become fruit sink limited at very low crop loads and extra vine photosynthates are diverted to vegetative structures. Conservative balanced pruning (33 + 11) consistently undercropped the vines (mean Y:PW = 7 and quartile range = 4.6 to 8.6) by limiting yields, even at large vine sizes. This is a possible pruning strategy for building vine size and improving the long-term productivity potential of young or stressed vineyards, but it is not a viable option for sustainable economic production. It is also consistent with other cool climate crop load studies which describe excessive canopy growth and the need for canopy management at Y:PW below 5. 17The conservative to moderate balanced pruning formulas used by Shaulis[START_REF] Shaulis | The interaction of resistant rootstock to the nitrogen, weed control, pruning, and thinning on the productivity of Concord grapevines[END_REF] promoted Concord crop loads that increased vine size each year. Since both vegetative and reproductive growth increased with vine size, the annual increase in vine size led to excessively large vines with increased leaf area, canopy density, and fruit zone shading. This led to the development and success of the Geneva Double Curtain divided canopy training and shoot positioning system to expose excessive leaf area to sunlight[START_REF] Shaulis | Response of concord grapes to light, exposure and geneva double curtain training[END_REF]). However, the current economic pressures in the Concord industry has, in general, led to chronic overcropping and a reduction in vine size below the point where canopy division is beneficial.Crop loads above 7.2 had a negative linear relationship with harvest juice soluble solids and seasonal change in vine size. In most other viticulture situations, where high fruit maturation is desired, this would be considered overcropping with incremental increases in Y:PW leading to greater reductions in juice soluble solids or later harvest dates and decreases in vine size and return crop potential. Liberal balanced pruning (66 + 11) had a mean Y:PW of 11 (quartile range = 8.3

  2. The seasonal change in vine pruning weight also increased with decreasing crop load but did not plateau. The crop load response to juice soluble solids and vine size was the same whether the crop load was adjusted through changes in pruning severity or by light pruning followed by mid-season fruit 21 thinning. There was a positive linear relationship between seasonal GDD and the Y:PW needed to reach 16 Brix, as well as between seasonal precipitation and the Y:PW required to observe no seasonal pruning weight change. The results of this study were used to develop a crop load management model for NY Concord which integrates the juice grape industry parameters for 16 Brix fruit, sustainable vine size and production, spring frost mitigation, seasonal environmental variation, and possible climate change impacts. Annually, the CLEREL research group tracks weather, Concord phenology, and Concord fresh berry weight on a set of sentinel vines and

Figure 1 Figure 2

 12 Figure 1 The effect of Concord pruning weight (PW) on fresh fruit yield at two variable node and two fixed node pruning severities. Vines were balanced pruned to retain 33 (open circles) or 66 (closed circles) nodes for the first 500g pruning weight plus an additional 11 nodes for each additional 500g or pruned to retain 100 (open triangles) or 120 (closed triangles) nodes per vine, regardless of pruning weight. Individual vine observations (n = 36 per year) from 2001-2004 were combined by pruning treatment and binned by 0.1 kg pruning weight (bars = ±SE). Fitted curves for 33 + 11: Yield = 11.8 -((11.8) * EXP(-1.1 * PW)), R 2 = 0.91; 66 +11: Yield = 14.9 -((14.9) * EXP(-1.3 * PW)), R 2 = 0.94; 100 nodes: Yield = 13.6 -((13.6) * EXP(-3.1 * PW)), R 2 = 0.87; 120 nodes:) * EXP(-3.6 * PW)), R 2 = 0.88.

Figure 3

 3 Figure 3 The effect of Concord Yield:Pruning weight (Y:PW) on the seasonal change in pruning weight from 2001-2004. Vines were grouped by pruning severity (33 + 11, 66 + 11, or 100 nodes: open symbols) or by 120 node pruning with mid-season fruit thinning (closed symbols). There were 108 pruning observations per year (36 per pruning level) and 75 thinning observations per year. Data was binned by crop load in intervals of 2 Y:PW. (bars = ±SE).

Figure 4

 4 Figure 4The general 4-year Concord crop load response on harvest juice soluble solids (A) and seasonal change in vine pruning weight (B). Pruning severity (open circles) or mid-season fruit thinning (closed circles) was used to generate a range of Y:PW. Each point is the mean of four years (n = 4). Pruning and thinning data were combined for the mean response. A broken stick model was fitted to illustrate the inflection point on the juice soluble solids response (see Equation1). For the pruning weight response, a linear model was used (see Equation2). Dashed lines indicate 95 % mean confidence limits.

Figure 5 Figure 7

 57 Figure 5 The effect of Concord Yield:Pruning weight on harvest juice soluble solids (A) and the seasonal change in vine pruning weight (B) from 1999 to 2009. Dashed lines are individual year trend lines and the solid lines are the mean response for all years. (n = 75 for 1999-2004, 2006, and 2010; 18 for 2005; and 42 for 2007-2009). The solid lines are the mean response (n = 11 years).

Table 1

 1 Weather and phenology from Cornell's Taschenberg Vineyard Laboratory, Fredonia, NY, from 1999-2009.

	2009. Mechanical cane pruning and crop adjustment decreases labor costs and
	maintaines fruit quality in New York "Concord" grape production. Horttechnology 19:247-
	253.
	Bates T, Dresser J, Eckstrom R, Badr G, Betts T, Taylor J. 2018. Variable-rate mechanical crop
	adjustment for crop load balance in "Concord" vineyards. 2018 IoT Vert Top Summit Agric
	-Tuscany, IOT Tuscany 2018:1-4.
	Bates TR. 2008. Pruning level affects growth and yield of New York Concord on two training
	systems. Am J Enol Vitic 59:276-286.

Bates TR. 2017. Mechanical crop control in New York "Concord" vineyards target desirable crop load levels. In Acta Horticulturae. pp. 259-264. International Society for Horticultural Science (ISHS), Leuven, Belgium.

Table 2

 2 Multivariate analysis of Yield:Pruning weight (Y:PW) and Treatment (pruning or thinning) on harvest juice soluble solids or the seasonal change in pruning weight in New York Concord from 2001-2004. Model used standard least squares and effect tests (Prob > F) are shown.

			Year
		2001	2002	2003	2004	All Years
	Juice soluble solids
	Y:PW	* a	***	***	**	***
	Experiment	***	***	*	NS	NS
	Y:PW * Experiment	*	NS	***	NS	NS
	Seasonal PW change
	Y:PW	*	***	***	***	***
	Experiment	NS b	**	NS	*	NS
	Y:PW * Experiment	NS	NS	NS	NS	NS

a *, **, *** indicate significance at p <0.01, <0.001, and <0.0001, respectively. b NS = not significant (p >0.01).

Table 3

 3 Multivariate regression effect tests of seasonal growing degree days or precipitation on the Yield:Pruning weight (Y:PW) needed to reach harvest juice soluble solids of 16 Brix or no seasonal change in vine pruning weight in New York grown Concord from 1999-2009.

			Y:PW for 16 Brix	Y:PW for no PW change
		F Ratio	Prob > F	F Ratio	Prob > F
	GDD	13.8	0.0099	0.1591	0.7019
	Precipitation	2.4	0.1700	10.66	0.0138
	GDD * Precip	0.04	0.8574	2.63	0.1483

This research was supported by the Lake Erie Grape Processors, the New York Wine and Grape Foundation, and the USDA-NIFA Specialty Crop Research Initiative Award No. 2015-51181-24393. The authors thank the late Dr. Nelson Shaulis for establishing the pruning level research at the Taschenberg Laboratory in Fredonia, NY, and acknowledge the effort of laboratory staff members Ted Taft, Paula Joy, Kelly Link, and Madonna Martin in vineyard maintenance and data collection.

26

Harvest delays and reduced crop potential for the following season are expected; however, vines will maintain balance in warmer and wetter than average seasons. This crop load recommended if mid-season fruit thinning is part of the management strategy. In cool and average seasons, the crop can be moderately thinned to maintain balance. In warm seasons, no thinning would be necessary.

> 25 Severely overcropped <15.7 JSS well below the 16 Brix standard and, if left unthinned, will still require a significant period of ripening after harvest has started. Vine size will be reduced by > 0.1 kg/vine (0.25 lbs/vine) with a lower future yield potential and a lower return crop. It requires excessive fruit thinning to achieve vine balance mid-season, which has been shown to cause canopy damage in Concord and negates the positive effects of fruit thinning on vine size/health. This level of crop load stress is not recommended. a Predicted Brix in an average season at a standard harvest of 30 -40 days after veraison. The given ranges reflect this spread of time. b An average season = 1455-1723 GDD (+/-1 st. dev. from the 11 -year GDD mean). Cool season < 1455 GDD, Warm season > 1723 GDD.