The anisotropic Kerr nonlinear refractive index of the beta-barium borate (\beta-BaB2O4) nonlinear crystal

We study the anisotropic nature of the Kerr nonlinear response in a beta-barium borate (\beta-BaB2O4, BBO) nonlinear crystal. The focus is on determining the relevant $\chi^{(3)}$ cubic tensor components that affect interaction of type I cascaded second-harmonic generation. Various experiments in the literature are analyzed and we correct the data from some of the experiments for contributions from cascading as well as for updated material parameters. We find that the Kerr nonlinear tensor component responsible for self-phase modulation in cascading is considerably larger than what has been used to date. We evaluate the impact of using such a cubic anisotropic response in ultrafast cascading experiments.


Introduction
Cascaded second-harmonic generation (SHG) describes the case where the frequencyconverted second harmonic (SH) is strongly phase mismatched, and this may create a Kerr-like nonlinear phase shift on the fundamental wave (FW). Since the early discovery of this cascaded self-action nonlinearity [1,2], cascaded SHG was for a long time overlooked until DeSalvo et al. measured a negative nonlinear phase shift from phase-mismatched SHG using the Z-scan method [3] (for an early comprehensive review on cascading, see [4]). The exciting promise of cascading nonlinearities is seen directly from the scaling properties of the Kerr-like nonlinear index it generates: n I 2,casc ∝ −d 2 eff /∆k, where d eff is the effective nonlinearity of the quadratic interaction. By simply changing the phase mismatch parameter ∆k we therefore have a Kerrlike nonlinearity that can be tuned in sign and strength. A particularly interesting and attractive property is that a self-defocusing Kerr-like nonlinearity n I 2,casc < 0 is accessible by having a positive phase mismatch ∆k > 0.
A popular crystal for cascading experiments is beta-barium-borate (β -BaB 2 O 4 , BBO), see e.g. [5][6][7][8][9][10][11][12][13][14][15]. It has a decent quadratic nonlinear coefficient, and because the crystal is anisotropic it can be birefringence phase-matched for type I (oo → e) SHG. For femtosecond experiments it has the important properties of a low dispersion, a high damage threshold, and a quite low Kerr self-focusing nonlinearity; the latter is important for cascading as we discuss later, because the material Kerr nonlinearity will compete with the induced cascading nonlinearity. The main goal of this paper is to analyze experiments in the literature where the Kerr nonlinear refractive index was measured in BBO [13][14][15][16][17][18][19][20][21], and for the fist time extract all four tensor components that are relevant to describe the anisotropic nature of the cubic nonlinearity.
While the anisotropy of the quadratic nonlinearity is extensively used for optimizing the phase-matching properties and maximizing the quadratic nonlinear coefficients in nonlinear crystals, the case is different when it comes to treating the cubic nonlinearities: the anisotropic nature has often been neglected when studying self-phase modulation (SPM) effects from the Kerr nonlinear refractive index, in contrast to the case of third-harmonic generation, as evidenced even in early studies [22,23]. Here we show which anisotropic nonlinear susceptibility components the experimental literature data represent. The experimental data will be analyzed and we correct the reported values for cascading contributions if necessary. We eventually obtain complete information of all four relevant tensor components for the BBO cubic nonlinear susceptibility. Our results show that the tensor component affecting the SPM of the ordinary wave, c 11 YYYY , is very well documented, and the corrected literature data agree extremely well in the near-IR with the popular two-band model [24]. We also perform our own experimental measurement of c 11 with a very accurate technique based on balancing out the self-focusing SPM from the Kerr effect with the self-defocusing effect from cascading (originally used in [15]). Both this measurement and the corrected data from the literature show that c 11 is substantially larger than what has been used so far in simulations in the literature, and an important consequence for cascaded SHG is that this reduces the range of phase-mismatch values that gives a total defocusing nonlinearity. In contrast experiments measuring the tensor components that affect the extraordinary wave interaction are more scarce, so new measurements are necessary in order to get more accurate and reliable values. We also assess the impact of using the anisotropic Kerr response in simulations of cascaded SHG, and it turns out that for BBO the optimal interaction angles are such that the anisotropy plays a minor role in the wave dynamics.
As to analyze the experiments properly and rule out any misunderstandings in the definitions, the foundation of the paper is a careful formulation of the propagation equations of the FW and SH waves under slowly-varying envelope approximation, see appendix A. These include the anisotropic nature of the frequency conversion crystal in formulating the effective cubic nonlinear coefficients. This theoretical background is important in order to analyze and understand the experiments from the literature that report Kerr nonlinearities for BBO. The analysis presented here should help understanding what exactly has been measured, and put the results into the context of cascaded quadratic soliton compression.

Background for cascading nonlinearities
In order to understand the importance of knowing accurately the Kerr nonlinear index, let us explain first how cascading works. The cascaded Kerr-like nonlinearity can intuitively be understood from investigating the cascade of frequency-conversion steps that occur in a strongly phase-mismatched medium [4]: After two coherence lengths 2π/|∆k| back-conversion to the FW is complete, and as ∆k = 0 the SH has a different phase velocity than the FW. Thus, the back-converted FW photons has a different phase than the unconverted FW photons, simply because in the brief passage when they traveled as converted SH photons the nonzero phase mismatch implies that their phase velocities are different. For a strong phase mismatch ∆kL ≫ 2π, this process repeats many times during the nonlinear interaction length, and in this limit the FW therefore effectively experiences a nonlinear phase shift parameterized by a higher-order (cubic) nonlinear index n I 2,casc given as [3] n I 2,casc = − where ω 1 is the FW frequency, n 1 and n 2 the linear refractive indices of the FW and SH, ∆k = k 2 − 2k 1 , and k j = n j ω j /c are the wave vectors. The cascaded cubic nonlinearity will compete with the material cubic nonlinearity (the Kerr nonlinear index n I 2,Kerr ), and the total nonlinear refractive index change experienced by the FW can essentially be described as ∆n = n I 2,tot I 1 = (n I 2,casc + n I 2,Kerr )I 1 (2) where I 1 is the FW intensity. Due to the competing material Kerr nonlinearity a crucial requirement for a total defocusing nonlinearity n I 2,tot < 0 is that the phase mismatch is low enough so that |n I 2,casc | > n I 2,Kerr (the so-called Kerr limit). When this happens, cascaded SHG can generate strong spectral broadening through SPM on femtosecond pulses, which can be compensated in a dispersive element. In a pioneering experiment by Wise's group femtosecond pulse compression was achieved this way [5]. Since the nonlinearity is self-defocusing the pulse energy is in theory unlimited as self-focusing effects are avoided [5], and it can even be used to heal small-scale and whole-scale self-focusing effects [6]. Another very attractive feature of the negative self-defocusing nonlinearity is that temporal solitons can be excited in presence of normal (positive) group-velocity dispersion, which means anywhere in the visible and near-IR. With this approach Ashihara et al. [7] achieved soliton compression of longer femtosecond pulses. Numerous experiments over the past decade have since been motivated by soliton compression of energetic pulses to few-cycle duration [9,11,[25][26][27]. Since the defocusing limit and the soliton interaction depend critically on the n I 2,Kerr value, knowing an exact value of the n I 2,Kerr coefficient is crucial.

Anisotropic quadratic and cubic nonlinearities in uniaxial crystals
In a uniaxial crystal, the isotropic base-plane is spanned by the crystal XY axes, and light polarized in this plane is ordinary (o-polarized) and has the linear refractive index n o . The optical axis (crystal Z-axis, also called the c-axis) lies perpendicular to this plane, and light polarized in this plane is extraordinary (e-polarized) and has the linear refractive index n e . In a negative uniaxial crystal like BBO, n o > n e . The propagation vector k in this crystal coordinate system has the angle θ from the Z-axis and the angle φ relative to the X-axis, cf. Fig. 1(a), and the e-polarized component will therefore experience the refractive index n e (θ ) = [n 2 o / cos 2 θ + n 2 e / sin 2 θ ] −1/2 , while the o-polarized light always has the same refractive index. Degenerate SHG can either be noncritical (type 0) interaction, where the FW and SH fields are polarized along the same direction, or critical (type I) interaction, where the FW and SH are cross-polarized along arbitrary directions in the (θ , φ ) parameter space. In type 0 the FW and SH are usually polarized along the crystal axes (θ = 0 or π/2) since this turns out to maximize the nonlinearity, and it is called noncritical because the interaction does not depend critically on the propagation angles (both nonlinearity and phase matching parameters vary little with angle).  1. (a) The definition of the crystal coordinate system XY Z relative to the beam propagation direction k. (b) Top view of the optimal crystal cut for type I oo → e SHG in BBO, which has φ = −π/2 and θ = θ c for perpendicular incidence of an o-polarized FW beam, and the e-polarized SH is generated through type I oo → e SHG. The specific value of the cut angle θ c depends on the wavelength and the desired application. Angle-tuning the crystal in the paper plane will change the interaction angle θ . Capital letters XY Z are traditionally used to distinguish the crystal coordinate system from the beam coordinate system xyz that has its origin in the k-vector propagation direction.
A great advantage of this interaction is that there is little or no spatial walk off (a consequence of the choice θ = 0 or π/2). In type I θ and φ values can often be found where phase matching is achieved. This interaction is very angle-sensitive, which is why it is called critical interaction.
The slowly-varying envelope equations (SVEA) for degenerate SHG are derived in Appendix A. There the quadratic and cubic "effective" anisotropic nonlinear coefficients were symbolically introduced. Below we show the expressions for these coefficients relevant for SHG in BBO, which belongs to the crystal class 3m, and in all cases we consider the reduced numbers of tensor components that come from adopting Kleinman symmetry (which assumes dispersionless nonlinear susceptibilities [28, Ch. 1.5]); thus BBO has 3 independent χ (2) tensor components and 4 independent χ (3) tensor components. An excellent overview of the quadratic and cubic nonlinear coefficients in other anisotropic nonlinear crystal classes is found in [21].
For an arbitrary input (i.e. FW) polarization various SHG processes can come into play (oo → o, oo → e, oe → e, oe → o, ee → e, and ee → o). This includes both type 0, type I and type II (nondegenerate SHG, where the FW photons are cross-polarized). Their respective d eff -values for a crystal in the 3m point group (which also includes lithium niobate) are [21] d eee eff = d 22 cos 3 (θ + ρ) sin 3φ + 3d 31 sin(θ + ρ) cos 2 (θ + ρ) + d 33 sin 3 (θ + ρ) (6) We remark here that the correct evaluation of the effective nonlinearity requires to take into account the spatial walk-off angle ρ = arctan[tan(θ )n 2 o /n 2 e ] − θ (for a negative uniaxial crystal) [29]. In BBO ρ is around 3 − 4 • for the typical angles used, which gives changes in d eff of a few percent, and as this is well within the standard errors on the experimental d i j values it is typically ignored for qualitative simulations. However, when we later evaluate the cascading contributions and compare them to the competing material Kerr nonlinearity, an exclusion of ρ would lead to a systematic error so we choose to include it in the analysis in Sec. 4.
BBO is usually pumped with o-polarized light as the quadratic nonlinearity is largest for this configuration (d 22 is the largest tensor component, see App. C), and because the oo → e interaction can be phase matched by a suitable angle θ . For such an interaction, the crystal is cut so φ = −π/2 as to optimize the nonlinearity (it is here relevant to mention that d 22 and d 31 has opposite signs in BBO). This has the direct consequence that oo → o interaction is zero. For this reason, when studying BBO pumped with o-polarized light in a crystal with φ = −π/2, Eq. (4) shows that d eff = d 31 sin(θ + ρ) − d 22 cos(θ + ρ).
Note that BBO has historically been misplaced in the point group 3, and in addition there has been some confusion about the assignment of the crystal axes (where the mirror plane of the crystal was taken parallel instead of perpendicular to the crystal X-axis) [29]. Unfortunately this means that even today crystal company web sites operate with d 11 as being the largest tensor component (in the 3m point group d 11 = 0) and supply crystals apparently cut with φ = 0 (because using the point group 3 combined with a nonstandard crystal axes definition means that d 11 cos 3φ must be maximized, while with the correct point group, 3m, and correct crystal axes assignments, d 22 is the largest tensor component and sin 3φ = 1 maximizes the effective nonlinearity). In order to sort out any confusion, the optimal crystal cut is shown in Fig. 1(b), and we hereby urge crystal companies to follow the standard notation.
If the Kerr nonlinearity is considered isotropic, we have χ      A mistake often seen in the early literature on type I interaction is to take χ XPM , and this mistake gives a 3 times too large XPM term. As shown above the identity χ only holds in the type 0 configuration, where the FW and SH have identical polarizations, but importantly it is not an identity that is restricted to isotropic nonlinearities as it holds for an anisotropic medium as well. Thus, the error made in the past for type I could either come from using directly the propagation equations for an isotropic medium and where two pulses with the same polarization interact, but it could also come from generalizing type 0 SHG propagation equations to type I, and forgetting the XPM properties for cross-polarized interaction.

Experimental conditions for the literature measurements
Most experiments aiming to measure the cubic nonlinearities have used the Z-scan method [38], which uses a Gaussian laser beam in a tight-focus limiting geometry to measure the Kerr nonlinear refractive index. The Z-scan method measures the transmittance of a nonlinear medium passing through a finite aperture placed in the far field as a function of the sample position (z) measured with respect to the focal plane. As the sample is Z-scanned (i.e. translated) through the focus of the beam, the lens effect from the Kerr nonlinear index change will change the amount of light recorded by the detector; this gives information about the intensity-induced nonlinear index change ∆n and thus the n I 2 coefficient, defined phenomenologically from the general expansion The intensity is rarely high enough to allow for any contributions but the n I 2 term. The superscript I underlines that the n I 2 parameter is the intensity-dependent nonlinear index, as one can also define an electric-field dependent index, typically as ∆n = 1 2 n E 2 |E| 2 . In the other sections of this paper we use n I 2,Kerr to denote the material Kerr nonlinear refractive index and n I 2,casc the Kerr-like cascading nonlinear index from cascaded SHG; this notation should suffice to avoid confusion with the linear refractive index of the SH, n 2 . As the result of a closed-aperture measurement can be influenced by contributions from multi-photon absorption, the aperture can be removed (open aperture scan). The n I 2 contributions then vanish, and only multi-photon absorption effects remain. Thus it is possible to separate the two contributions.
There are some issues with the Z-scan method that may affect the measured nonlinear refractive index. If the repetition rate is too high, there are contributions to the measured n I 2 from thermal effects as well as two-photon excited free carriers [39] (for more on these issues, see e.g. [40]). Similarly, a long pulse duration can also lead to more contributions to the measured n I 2,Kerr than just the electronic response that we aim to model, and in particular the static (DC) Raman contribution is measured as well. However, for BBO the fraction of the delayed Raman effects is believed to be quite small, so we will assume that it can be neglected making the measured nonlinearities correspond to the electronic response.
In the following we will analyze measurements in the literature in order to extract the four BBO cubic tensor components in Eqs. (7)- (9). When converting from the nonlinear susceptibilities c µm to the intensity nonlinear refractive index we can use Eq. (53). There is a small caveat here, because the nonlinear index is defined for SPM and XPM terms, where only two waves interact through their intensities, while the c µm coefficients are generally defined for four-wave mixing. Thus, it does not always make sense to define a nonlinear refractive index (consider e.g. the c 10 component) until the total effective susceptibility is calculated, e.g. through Eqs. (7)- (11). Thus, it is safer to keep the susceptibility notation as long as possible when calculating the effective nonlinearities. The exception is when measuring the nonlinear refractive index of the oooo SPM interaction because the relevant nonlinear tensor component that is measured is c 11 , see Eq. (7), and thus the connection between c 11 (i.e. χ (3) oooo ) and n I 2,Kerr is unambiguous. In the experiments in the literature the cascaded quadratic contributions were often forgotten or assumed negligible. The exact relations for the evaluating cascading nonlinearities are derived in Appendix B. We will address these issues below by in each case calculating where possible the cascading nonlinearity.

The c 11 tensor component
In this section we present measurements using o-polarized input light, which means that the c 11 tensor component is accessed by measuring the SPM component of the o-polarized light χ To our knowledge, the first measurement of any of the cubic susceptibilities was performed by Tan et al. [13]. They investigated the cascaded nonlinearity in a 10 mm BBO crystal cut at θ = 22.8 • for phase matching at 1064 nm (they claim also that the crystal is cut with φ = 0, which is probably due to the confusion related to placing BBO in the point group 3 instead of 3m, cf. the discussion on p. 7). Pumping with 30 ps o-polarized pulses from a Nd:YAG modelocked laser at two nearly degenerate wavelengths (a strong wave at λ a = 1.064 µm and a weaker wave at λ b = 1.090 µm), they rotated the BBO crystal around ∆k = 0 in a type I oo → e SHG setup, where both the SH 2ω a of the strong beam was generated and four-wave mixing of the two pumps ω c = 2ω a − ω b was generated. By recording the impinging and generated pulse energies vs. angle they found by a fitting analysis the coefficient g = 220, where g =χ was the only free parameter in the fit. For θ = 22.8 • we have ρ = 3.2 • and d eff = 1.99 pm/V, so this gives the measured nonlinear susceptibilityχ They also performed a Z-scan measurement, and the setup was calibrated towards a BK7 glass sample. They foundχ BK7 ≃ 1.4. By using the value they proposeχ We report these results with a bar because their definition of the cubic susceptibility is different than what we use: they define the SPM polarization as P int E |E | 2 /2, which replaces our definition Eq. (47), so the connection between them isχ int = 6.08 × 10 −22 m 2 /V 2 . In our notation their cascading measurements gave (14) corresponding to n I 2,Kerr (ω 1 ; ω 1 ) = 6.03 ± 0.52 · 10 −20 m 2 /W. Their Z-scan measurement gives corresponding to n I 2,Kerr (ω 1 ; ω 1 ) = 4.18 ± 0.54 · 10 −20 m 2 /W. The error bars are indicative as they do not report the errors in their measurements, but we assume 5% error on g and conservatively estimate an error of 10% on the Z-scan result. A question remains concerning the Z-scan result: they do not report the angle θ at which the measurement was done so we cannot asses the amount of cascading contributions they had in this measurement, but given the fact that the c 11 value is quite small it seems likely that cascading had a contribution.
Hache et al. [14] used a 1 mm BBO cut at θ c = 29.2 • with 800 nm 100 fs pulses from a MHzrepetition rate Ti:sapphire oscillator. The pump was o-polarized allowing for oo → e cascading interaction. They carried out Z-scan measurements both close to the phase-matching point θ c as well as "far from SHG phase matching". In the latter case they found which corresponds to c 11 = 4.36 ± 0.97 · 10 −22 m 2 /V 2 . They claim that cascading contributions did not contribute for this measurement, but whether this is true cannot be judged with the information at hand, as they did not specify the angle they used for this measurement. DeSalvo et al. [16] used the Z-scan method with a single-shot 30 ps pulse at 1064 nm, which had k parallel to the optical axis c, so θ = 0 and consequently ρ = 0 as well. Thus, (a) the pump is always o-polarized no matter how the input beam is polarized, and (b) the input polarization determines the angle φ . Since the specific orientation of the crystal with respect to the input polarization was not reported, we take it as unknown. The relevant Kerr contribution at the pump wavelength is oooo interaction, which is given by Eq. (7) as simply one tensor component c 11 . The experiment reported n E 2 = 11 ± 2 · 10 −14 esu at 1064 nm, where the nonlinear index change is defined as ∆n = 1 2 n E 2 |E| 2 , and using Eq. (C1) in [30] this corresponds to χ (3) = 2.70 ± 0.49 · 10 −22 m 2 /V 2 . The cascading contributions are from oo → e and oo → o processes. Since θ = 0 they have the same phase mismatch values as e-polarized light in this case has the same refractive index as o-polarized light. At 1064 nm the value is ∆k ooe = ∆k ooo = 234 mm −1 . Taking into account the two cascading channels Eqs. (4) and (3), and using Eq. (55), we get the total contribution χ i.e. independent on the propagation angle φ . For these calculations we used d 22 = −2.2 pm/V. The d eff values are typically reported with a 5% uncertainty [42], giving a 7% uncertainty on the cascading estimate. This means that when correcting for the cascading contributions we get c 11 = 4.65 ± 0.51 · 10 −22 m 2 /V 2 , λ = 1.064 µm (17) corresponding to n I 2,Kerr (ω 1 ; ω 1 ) = 4.80 ± 0.53 · 10 −20 m 2 /W. At 532 nm n E 2 = 21 ± 4 · 10 −14 esu was measured, corresponding to χ (3) = 5.22 ± 0.99 ·10 −20 m 2 /V 2 . The cascading is smaller here, using ∆k ooe = ∆k ooo = 1, 933 mm −1 we get χ  (18) corresponding to n I 2,Kerr (ω 1 ; ω 1 ) = 5.87 ± 1.00 · 10 −20 m 2 /W. At 355 nm n E 2 = 14 ± 3 · 10 −14 esu was measured, corresponding to χ (3) = 3.54 ± 0.67 · 10 −22 m 2 /V 2 . The cascading is small, using ∆k ooe = ∆k ooo = 11, 736 mm −1 we get χ  (19) corresponding to n I 2,Kerr (ω 1 ; ω 1 ) = 3.81 ± 0.74 · 10 −20 m 2 /W. We must mention that the cascading value is quite uncertain: the SH wavelength is here 177.5 nm, which is very close to the UV poles in the Sellmeier equations (the fit in [43] puts the poles at 135 nm and 129 nm for oand e-polarized light, respectively). Finally, at 266 nm they measured n E 2 = 1 ± 0.3 · 10 −14 esu.
Moses et al. [15] pumped a BBO crystal at 800 nm with 110 fs o-polarized light from a 1 kHz Ti:Sapphire regenerative amplifier. By angle-tuning the crystal they achieved zero nonlinear refraction: this happens when the cascading from the oo → e interaction exactly balances the Kerr nonlinear refraction. Investigating a range of intensities they determined this critical phase mismatch value to ∆k c = 31 ± 5 π/mm. Using this value and assuming that it corresponds to n I 2,tot = 0 then a reverse calculation though Eq. (1) gives n I 2,Kerr = −n I 2,casc . Using the value d eff = 1.8 pm/V they inferred the value n I 2,Kerr = 4.6 ± 0.9 · 10 −20 m 2 /W (the main sources of the uncertainty are to determine the precise phase-mismatch value as well as the uncertainty on d eff ). Since o-polarized pump light was used, this contains only the c 11 tensor contribution from an oooo interaction. However, the value d eff = 1.8 pm/V they used to infer this nonlinearity was probably taken too low [45]. Let us therefore estimate it again using a more accurate value: ∆k = 31 π/mm is achieved at 800 nm with θ = 26.0 • , giving ρ = 3.6 • and d eff = 2.05 pm/V. With this corrected value we get n I 2,Kerr (ω 1 ; ω 1 ) = 5.54 ± 0.98 · 10 −20 m 2 /W, λ = 0.8 µm (28) which corresponds to c 11 = 5.41 ± 0.95 · 10 −22 m 2 /V 2 .
In the same spirit, we performed an experiment using the technique employed by Moses et al. The crystal was a 25 mm BBO from Castech with a 10x7 mm 2 aperture, cut so θ = 21 • and φ = −π/2, and antireflection coated. The crystal was pumped with pulses from a 1 kHz Ti:Sapphire regenerative amplifier followed by a commercial OPA, giving λ = 1.032 µm pulses with around 40 µJ pulse energy. The pulses were Gaussian with 39.3 nm FWHM bandwidth and an intensity autocorrelation trace revealed a Gaussian pulse with an intensity FWHM of 46 fs; the pulses were therefore nearly transform-limited Gaussian pulses (the time-bandwidth product was 16% above transform limit). The pump was collimated with a telescope before the crystal and the spot size was measured to 0.5 mm FWHM, and the intensity was controlled by a neutral-density filter. We chose to pump at this wavelength because the GVM-induced selfsteepening is smaller than at 800 nm, and besides the transition from negative to positive effective nonlinearity occurs in the so-called stationary regime, where the cascading nonlinearity is non-resonant (see also Sec. 4.5 later). We first located the phase-matching angle and then angletuned the crystal until any visible nonlinear components of the spectrum had vanished (i.e. it returns close to its input state) as n I 2,casc + n I 2,Kerr ≃ 0 in this range. By fine-tuning (with 1/6 • precision) θ in this range we measured the trends shown in Fig. 2(a): this plots the bandwidth 10 dB below the peak value of the recorded spectra for two different intensities, and the curves are seen to cross each other at ∆k ≃ 88 mm −1 . We investigated numerically whether this is a good measure of the transition from negative to positive total nonlinearity; for the parameters used in the numerics the critical point was calculated analytically to be ∆k c = 88.5 mm −1 (neglecting XPM effects, see later). Fig. 2(b) shows the equivalent numerical results, and clearly the two different intensities cross each other very close to the theoretical transition. Two comments are in order: firstly, the model includes XPM effects that give a self-focusing contribution [37] (see also discussion later in Appendix B, where the cascading quintic term is calculated) that for the chosen intensities is quite minimal. We can judge the XPM impact by turning off XPM in the  code: the transition now occurs practically at the theoretical value, which is a bit higher than the result with XPM. Thus, we can expect that XPM gives a small correction on the order of −1 mm −1 . Secondly, the reason why the spectrum does not return to the input value is that when cascading cancels Kerr SPM effects there are still self-steepening effects that remain (see the detailed experiments by Moses et al. [15], and also the recent discussion by Guo et al. [46]). Plots (c) and (d) show the spectra right at the transition, and the red edge is quite similar for the two intensities while the blue edge is always higher for the larger intensity. We should mention that we found similar transition "indicators" at very similar ∆k values by observing the -20 dB bandwidth, FWHM bandwidth, in tracking the wavelength position of the spectral edges, in tracking the center wavelength etc. Especially the latter is the most precise indicator   [15]. The theoretically predicted electronic nonlinearity is calculated with the 2-band model [24]. The average value curve was calculated through a weighted mean of the Miller's delta from all data, except the UV measurements below 400 nm, and the shaded areas denoted "σ " and "2σ " represent one and two standard deviations, respectively. numerically as cascading gives a very pronounced blue-shift of the spectrum -see graph (e)but experimentally the peak was too erratic in the fine-tuning around the transition, which is probably related to the fact that the pump is generated as the SH of the idler of our OPA. We conclude that after correcting for a small XPM contribution the critical transition is found the be ∆k c = 89 ± 5 mm −1 , where the uncertainty lies mostly in finding a suitable transition value to track between the different intensities and phase mismatch values. This critical value corresponds to θ = 18.5 • , ρ = 2.7 • and d eff = 2.08 pm/V, so calculating the cascading nonlinearity and employing the n I 2,casc + n I 2,Kerr = 0 hypothesis, we get n I 2,Kerr (ω 1 ; ω 1 ) = 4.87 ± 0.44 · 10 −20 m 2 /W, λ = 1.032 µm (29) which corresponds to c 11 = 4.73 ± 0.43 ·10 −22 m 2 /V 2 . The uncertainty might be reduced with a longer pump pulse as self-steepening effects are reduced and it therefore will be easier to compare the nearly-linear spectrum at ∆k c (where SPM is nulled and only XPM and self-steepening effects contribute) to the input spectrum, but a drawback in using a longer pulse could lie in the challenge of observing small changes in a narrower spectrum. All in all, we must emphasize that this method is quite accurate for determining the Kerr nonlinearity, especially compared to the 20-30% uncertainty often seen with other techniques. This is partially because we do not rely on any absolute pump intensity, which can be a large uncertainty, and because the quadratic nonlinearities and the linear dispersion properties of BBO both are very well known. This completes the c 11 measurements. A summary of the data is shown in Fig. 3, together with the predicted electronic nonlinearity from the two-band model (2BM) [24]. The experi-mental near-IR data are amazingly well predicted on an absolute scale by the 2BM; we must here clarify that for the 2BM we chose the material constant K = 3100 eV 3/2 cm/GW, which was found appropriate as a single material parameter for dielectrics [16], and by using the BBO band-gap value E g = 6.2 eV [16]. Thus, in practice there are no free parameters in the model, which underlines the incredible agreement obtained on an absolute scale (and not just the translation from one wavelength to another). We mention here that in DeSalvo et al. [16] they proposed to "rescale" the bandgap do obtain a better fit with the experimental data (mainly because the two-photon absorption values β were not accurately predicted with 6.2 eV), and eventually suggested using 6.8 eV instead of 6.2 eV for BBO. This gave a much better β agreement at 266 nm, while the 352 nm value was still off. Here we see that the corrected DeSalvo n I 2,Kerr near-IR values actually agree extremely well with the 2BM n I 2,Kerr value when using E g = 6.2 eV, while the agreement is less accurate with E g = 6.8 eV. A common issue whether one uses 6.2 or 6.8 eV is that the predicted strong enhancement just below 400 nm due to two-photon absorption of the FW is not reflected in the experimental data. This is of little importance for our purpose, because in cascaded SHG the FW wavelength is usually never below 600 nm. In order to extract a suitable value for the c 11 coefficient, we calculated first the Miller's delta from all the data. It is defined as [35, 36] where χ j (ω 0 ) = n 2 j (ω 0 ) − 1 is the linear susceptibility. We here excluded the UV measurements below 400 nm as they obviously are not represented well by Miller's scaling. Since Miller's delta should be independent of wavelength this should allow us to calculate an average value based on the measurements at different wavelengths. We confirmed this hypothesis by checking that that Miller's delta vs. measurement wavelength could be fitted to a line with a near-zero slope, before we calculated the weighted average (with a weight given by σ −2 j where σ j is the estimated error of the jth measurement) as ∆ 1111 = 52.3 ± 7.7 × 10 −24 m 2 /V 2 (31) where the standard deviation was calculated using an unbiased weighted average. This value corresponds to c 11 = 5.0 · 10 −22 m 2 /V 2 and n I 2,Kerr (ω 1 ; ω 1 ) = 5.1 · 10 −20 m 2 /W at 800 nm. We see in Fig. 3 that this average is a quite good representative of the experimental data; all except two are within two standard deviations. We also note that the value we measured in this work fits extremely well with this proposed average value.

Other tensor components
Banks et al. [21] measured the c 10 tensor component of BBO using third-harmonic generation (THG) measurements at 1053 nm using 350 fs pulses from a 10 Hz Ti:Sapphire regenerative amplifier. The propagation angle was θ = 37.7 • , and from the data recorded while φ was varied around zero they extracted C 10 = −6 · 10 −24 m 2 /V 2 using the quadratic nonlinear coefficients of Shoji et al.
The uncertainty of the measurement was not reported, but for a similar fit it was reported to be 15%, which is what we used above. They also measured a three times larger value by using different quadratic nonlinearities, namely d 15 = 0.16 pm/V (as opposed to d 15 = 0.03 pm/V). Banks et al. [21] also measured the mixture 1 3 C 11 cos 2 θ m + C 16 sin 2 θ m = 4.0 ± 0.2 · 10 −23 m 2 /V 2 at θ m = 47.7 • . The TH walk-off angle is ρ = 4.5 • , and replacing in the above expression θ → θ + ρ will not change significantly the following numbers; therefore we choose to ignore the walk-off contribution because the fit done by Banks et al. also ignored it. If we now use Eq. (31), and the apply Miller's scaling to get the c 11 component appropriate for ω + ω + ω → 3ω interaction, see Eq. (30), we get c 11 (−3ω; ω, ω, ω) = c 11 (THG@λ = 1.053 µm) = 5.24 · 10 −22 m 2 /V 2 , and therefore c 16 (THG) = 1.47 ± 0.34 · 10 −22 m 2 /V 2 , λ = 1.053 µm (33) A small caveat must be noted: what they measured was related to the χ (3) (−3ω; ω, ω, ω) coefficients, because they investigated the yield of the third harmonic with respect to the pump, and decoupled the cascade yield (from multistep SHG mixing, i.e. ω + ω → 2ω followed by ω + 2ω → 3ω) to get the pure cubic nonlinear contribution. However, it is not sure that χ (3) (−3ω; ω, ω, ω) = χ (3) (−ω; ω, −ω, ω), i.e. that the THG nonlinearity is the same as the SPM nonlinearity; this was discussed, e.g., in [47]. Therefore using these values to model the SPM and XPM effects in BBO is an approximation. Sheik-Bahae and Ebrahimzadeh [17] used the Kerr-lens autocorrelation method [48] to measure n I 2,Kerr = 3.65 ± 0.6 · 10 −20 m 2 /W (χ (3) = 3.43 ± 0.56 · 10 −22 m 2 /V 2 ) using a 76 MHz Ti:Sapphire oscillator giving 120 fs pulses at λ = 0.850 µm. The BBO crystal was cut for phase matching at 800 nm, i.e. with θ = 29.2 • , and consequently ρ = 3.9 • . The crystal cut was optimized for SHG, so φ = −π/2. The pump pulse was e-polarized as to ensure lack of phase matching and thereby no cascading. We will now assess whether this assumption is fulfilled. The cascading contributions that might occur are ee → o and ee → e, however the former has zero nonlinearity when φ = −π/2, cf. Eq. (5). The phase mismatch for the latter interaction is ∆k eee = 396 mm −1 for θ = 29.2 • . The cascading contribution at φ = −π/2 then becomes χ Let us now understand what tensor components this value represents. By using e-polarized light several tensor components come into play for the cubic nonlinearity: the effective Kerr nonlinearity experienced by the pump is an eeee interaction, which is given by Eq. (11). For θ = 29.2 • the TH walk-off angle is ρ = 4.3 • , and with φ = −π/2 the effective nonlinearity is χ (3) eff (ω 1 ; ω 1 ) = −1.28c 10 + 0.48c 11 + 1.27c 16 + 0.093c 33 (34) Four tensor components appear, and we note that we here have a measurement that involves the c 33  An exciting consequence of this result is that it points towards a negative, self-defocusing Kerr nonlinearity. However, as indicated the particular value is very uncertain (the uncertainty was calculated with error propagation rules), which mainly stems from the low prefactor in front of  the c 33 term in Eq. (34). Apart from that, it also relies on three other separate measurements, making it very sensible to the c 11 , c 10 and c 16 values used. Here we mention that the large uncertainty of c 10 reported in [21] plays a role; with the various possibilities reported there for the c 10 value we always get a negative c 33 , and its value can vary by a factor of two from that reported in Eq. (35).

Summary of experiments
A summary of the data reported in the literature along with our corrected or updated values is shown in Table 1. This table also gives an overview of the crystal angles, pulse duration, repetition rates, wavelengths and the tensor components accessed in the measurements. Our analysis of the experiments [14][15][16][17][19][20][21] points towards using the Kerr nonlinearities summarized in Table 2. Only the susceptibilities are reported in order to minimize the errors when using the values for composite nonlinear indices. We have there also indicated the Miller's delta calculated from Eq. (30). Since the experiments are performed at different wavelengths, using the ∆ i jkl values makes it easier to evaluate a linear combination of the nonlinear coefficient at some particular wavelength. [Another popular model for scaling the cubic nonlinearity is the Boling-Glass-Owyoung (BGO) model [50], but we did not see a big difference in using Table 2. Proposed nonlinear susceptibilities for the BBO anisotropic Kerr nonlinearity. Note the c 16 value is deduced from the c 11 and c 10 coefficients, and the c 33 value is deduced from the c 11 , c 10  that model compared to Miller's rule. Besides, extending the BGO model to anisotropic nonlinearities is not straightforward, so we prefer to use Miller's rule for frequency scaling.] The values are reported with more significant digits than supported by the uncertainties, but this is done on purpose so it is easier to cross-check the coefficients as well as frequency scaling them.

Implications for cascaded pulse compression in BBO
In context of cascaded quadratic nonlinearities by far most important component is the FW SPM coefficient χ eff (ω 1 ; ω 1 ), and in a type I oo → e configuration it is given by the c 11 component. In some of our previous work related to BBO [37,51,52], we assumed an isotropic Kerr nonlinearity, and have used the value n I 2,Kerr (ω 1 ; ω 1 ) = 3.65 · 10 −20 m 2 /W taken from [17], and in a later publication, we have used the value given by Eq. (28). An implication of this larger value is that the "compression window" [37] becomes smaller. With this we imply a range of phase-mismatch values, where compression is optimal. We first need ∆k to be small enough so n I 2,tot < 0, i.e. |n I 2,casc | > n I 2,Kerr (ω 1 ; ω 1 ). This is the upper end of the window, also denoted the Kerr limit. Note that as mentioned in Sec. 3 the Kerr limit is quite sensitive to whether one includes the spatial walk-off angle ρ in calculating the d eff -value in the cascading expression. Once below the Kerr limit, decreasing ∆k strengthens the total defocusing nonlinearity. However, at a certain point the group-velocity mismatch (GVM) becomes too strong and the compression quality is strongly reduced. Several things happen: (a) the GVM-induced self-steepening term [9,46] is increased as it scales as d 12 /∆k, where d 12 is the GVM parameter, so the compressed pulse experiences a strong pulse-front shock. (b) When ∆k becomes too low (specified more accurately below) the SH will experience a resonant phase matching condition of a sideband frequency to the center frequency of the FW spectrum [53]. This is damaging to pulse compression because in the cascading process effectively this is the cascading bandwidth felt by the FW, and in the transition from the nonresonant to the resonant regime, the bandwidth essentially goes from being octave spanning to becoming resonant and thereby narrow [27]. The criterion for being in the nonresonant regime, also denoted the stationary regime, can in the simple case where only second-order dispersion is considered [51] be expressed as d 2 12 − 2k 2 (ω 2 ) is the group-velocity dispersion (GVD) coefficient of the SH. Thus, in the case where the SH GVD is normal (positive) we have the stationary (nonresonant) regime with broadband cascading when ∆k > ∆k sr ≡ d 2 12 /[2k (2) 2 (ω 2 )]. Below this threshold, in the so-called nonstationary regime, the cascaded nonlinearity is as mentioned resonant: the poor bandwidth implies that there is no possibility to achieve few-cycle duration [30]. At the same time the compressed pulse quality is low as there is a strong pulse shock front: this stems from the d 12 /∆k ratio being large, and this term controls the cascading- In order to excite solitons ∆k must be kept below the Kerr limit (red line). Optimal compression occurs when the cascaded nonlinearities dominate over GVM effects (∆k > ∆k sr , above the black line). Note that ∆k sr is calculated for the full dispersion case, and that for λ 1 > 1.49 µm the FW GVD becomes anomalous. We have also indicated the operation wavelengths of Cr:forsterite, Yb and Ti:Sapphire based amplifiers. The Kerr limit employs Miller's rule to calculate the nonlinear quadratic and cubic susceptibilities at other wavelengths, case (1) corresponds to the 'old' Kerr value n I 2,Kerr (ω 1 ; ω 1 ) = 3.65 · 10 −20 m 2 /W@850 nm (taken from [17]) and case (2) corresponds to the Kerr value proposed in this work ∆ 1111 = 52.3 m 2 /V 2 , corresponding to n I 2,Kerr (ω 1 ; ω 1 ) = 5.06 · 10 −20 m 2 /W@850 nm. (b) Numerical simulation of the case marked with 'x' in (a): a 50 fs@1030 nm I in = 500 GW/cm 2 pulse propagating in a 25 mm BBO crystal with θ = 19.1 • , ρ = 2.8 • and φ = −90 • (∆k = 80 mm −1 ). The simulations used the plane-wave SVEA equations (50)-(51) including full dispersion and extended to include self-steepening, and case (1) assumes an isotropic Kerr nonlinearity and n I 2,Kerr (ω 1 ; ω 1 ) = 3.65 · 10 −20 m 2 /W@850 nm, while (2) uses the anisotropic coefficients of Table 2, and Eq. (30) to calculate the nonlinear coefficients at 1030 nm.
induced self-steepening, as mentioned above. As a consequence one can only achieve a decent pulse quality by keeping very low soliton orders (roughly below 2) [30]. The implications of the larger nonlinear SPM Kerr nonlinearity for the FW is visualized in Fig. 4(a). With the lower value, n I 2,Kerr (ω 1 ; ω 1 ) = 3.65 · 10 −20 m 2 /W@850 nm, a large compression window is predicted, almost even encompassing Ti:Sapphire laser wavelengths. With the new larger value, ∆ 1111 = 52.3 m 2 /V 2 corresponding to n I 2,Kerr (ω 1 ; ω 1 ) = 5.06 · 10 −20 m 2 /W@850 nm, the compression window starts opening around 900 nm. To judge the impact of using this larger Kerr value, Fig. 4(b) shows a simulation at 1030 nm, where the phase mismatch is chosen to lie eff (ω 1 ; ω 1 ). The range indicated corresponds to the θ range in (a). The dashed lines indicate the isotropic limit. All nonlinear coefficients are scaled to other wavelengths using Miller's rule, and we used the SHG nonlinearities at 1064 nm from App. C as well as the cubic nonlinear parameters listed in Sec. 5. We also took φ = −π/2. inside the compression window using the new larger Kerr value. After 25 mm of propagation the short (50 fs) and intense (I in = 500 GW/cm 2 ) input pulse is compressed somewhat due to the formation of a self-defocusing soliton, and the spectrum shows SPM-like broadening and even formation of a soliton-induced Cherenkov wave (dispersive wave) around 2900 nm [12]. When instead using the 'old' lower Kerr value (dashed lines) the result is changed quite a lot: this is because |n 2,tot | is now much larger so the intensity gives a larger effective soliton order N 2 eff ∝ |n 2,tot |I in (the specific numbers are N eff = 3.4 in case (1) and N eff = 1.7 in case (2). The soliton is therefore compressed more and earlier in case (1), so after 25 mm soliton splitting has already occurred and in the spectrum the broadening is more pronounced and the dispersive wave is much stronger. These results clearly show how sensible cascading is to the value of the FW Kerr SPM parameter.
The experiments [5,7] carried out at 800 nm were performed well into the nonstationary region, which was necessary to achieve a defocusing nonlinearity. No few-cycle compressed solitons were observed in [7], as the intensity had to be kept low in order to avoid severe GVMinduced self-steepening. (Note that the simulations in [7] used χ (3) = 5 · 10 −22 m 2 /V 2 @800 nm, i.e. identical to the value we suggest.) Instead the experiment carried out at 1250 nm [10] was done in the stationary regime, and indeed a few-cycle soliton was observed.
A question is whether pumping with e-polarized light could give any advantages. As the ee → o interaction is heavily phase mismatched in a negative uniaxial crystal, the main cascading channel would be the ee → e interaction. As expected this interaction can never be phase matched. The maximum d eff is for θ = 0, and this turns out to give the maximum cascading strength as well. However, our calculations show that it is too small to overcome the n I 2,Kerr (ω 1 ; ω 1 ) nonlinearity, and the total effective cubic nonlinearity would be focusing. If θ = π/2 we get d eff = d 33 , which for BBO is very small (this is typical for borates, while the niobates instead have very large d 33 components). Nevertheless, for θ = π/2 the c 33 component determines the SPM coefficient, and as we found it could be negative. The total effective cubic nonlinearity would therefore be de-focusing. Investigating a BBO crystal pumped with θ = π/2 using e-polarized light would therefore be an interesting next step to resolve this issue.
We should finally mention that the error made in assuming an isotropic response for the cascaded soliton compression is not large: most importantly, the crucial parameter is the FW SPM coefficient, which as we saw is the same in the isotropic and in the anisotropic cases for this type I interaction. The XPM and SH SPM coefficients change when anisotropy is taken into account, but they play a minor role and rather tend to perturb the result more than shape it. In order to make a more quantitative statement, we calculated the θ -ranges where self-defocusing solitons can be excited in BBO. A necessary (but not sufficient [37]) criterion for this to happen is ∆k > 0 as well as ∆k low enough for |n I 2,casc | > n I 2,Kerr (ω 1 ; ω 1 ). The range is shown vs. FW wavelength in Fig. 5(a). Inside this θ range we can now for each wavelength calculate the anisotropic Kerr nonlinearities. As we know the FW SPM coefficient does not change with θ (but it does change across the wavelength range shown, which we here estimate using Miller's rule), so we can normalize the results to this value. We then get the results shown in Fig. 5(b). The XPM term is seen to lie close to the isotropic value 1/3 (found by taking χ In summary, assuming an isotropic response in BBO is a good approximation for the FW XPM term, mainly because the working θ range lies close to zero (the isotropic limit). This result cannot be generalized to other nonlinear crystals, where the working range might lie closer to π/2; an example is lithium niobate that operates quite close to θ = π/2 in the type I configuration, see e.g. [30]. For the SH SPM term its value is quite far from the isotropic case, but off all the parameters this is the least important one because the SH SPM is negligible in the cascading limit compared to FW XPM and in particular FW SPM. To confirm this we checked that the simulations of case (2) in Fig. 4(b) were almost identical results when using isotropic Kerr nonlinearities. A similar conclusion was drawn in our recent work, where the anisotropy was also taken into account [46].

Summary
We have analyzed the Kerr nonlinear index in BBO from a number of experiments [14][15][16][17][19][20][21], and we argued that in many of them contributions from cascaded SHG nonlinearities need to be taken into account. We also performed a very accurate measurement of the main tensor component c 11 by balancing the Kerr nonlinearity with a negative cascading nonlinearity. Our analysis points towards using the Kerr nonlinearities summarized in Table 2. They encompass all four relevant tensor components for describing the anisotropic nature of the BBO Kerr nonlinearity under Kleinman symmetry. The most reliable value is the c 11 component, relevant for describing the oooo SPM coefficient, while the other three are more uncertain: the c 16 and c 10 values are measured for THG, so using them for describing SPM and XPM effects is an approximation. Finally, the c 33 coefficient is very uncertain as it was measured in a non-ideal way (since the experiment used a crystal angle where its importance is very small) and on top of that all the three other coefficients in the table were used to deduce its value. In all cases the measured coefficients had contributions from both instantaneous electronic and delayed Raman effects, but we assumed that the Raman contributions to the measured nonlinearities were vanishing.
The proposed c 11 value was calculated as a weighted average over 14 different experiments in the 532-1064 nm range, and is therefore expressed as a wavelength-independent Miller's delta ∆ 1111 = 52.3 ± 7.7 × 10 −24 m 2 /V 2 , cf. Eq. (31); at 800 nm this value corresponds to c 11 = 5.0 × 10 −22 m 2 /V 2 and n I 2,Kerr = 5.1 · 10 −20 m 2 /W, which is larger than typical values used in the literature. We stress that the c 11 value we measured in this work agrees extremely well with this proposed value, and we believe that this method -where the focusing Kerr nonlinearity is nulled by an equivalent but defocusing cascading nonlinearity -is quite precise and can be used in many other contexts as an alternative to Z-scan measurements. When plotting the (corrected) Kerr nonlinear refractive indices vs. wavelength, see Fig. 3, we came to a surprisingly good absolute agreement with the 2-band model [24], except in the short end of the visible range.
Finally we showed the predicted consequences for cascaded femtosecond pulse compression exploiting type I interaction in BBO when using this larger value: the operation range where few-cycle pulse compression can be achieved will shift to longer wavelengths, roughly above 1.0 µm. We also showed that the anisotropic XPM coefficient of the FW is quite similar to the isotropic one, but this relies on the particular range of phase-mismatch rotation angles used in BBO, so this result cannot be generalized to other crystals exploiting type I interaction. Instead the anisotropic SH SPM coefficient is quite different from the one found assuming an isotropic Kerr nonlinearity, but since in cascading the SH is weak this difference will amount to much in a simulation of cascaded SHG.

Acknowledgments
The

A. The propagation equations under the slowly-varying envelope approximation
We have previously touched upon the issue of anisotropic cubic nonlinearities in quadratic nonlinear crystals [30], where we focused on type I interaction and the measurements in the literature for the lithium niobate crystal. Here we formulate the complete propagation equations that hold for both type 0 and type I, and discuss the anisotropic cubic tensor components relevant for the BBO crystal class.
Let us show briefly how to derive the basic slowly-varying envelope equations (SVEA) for degenerate SHG (i.e. where the FW photons are degenerate), thus describing both type 0 and type I SHG. We do this in order to avoid confusion about how to define the relation between the electric field and the susceptibility tensors. Working in mks units, from Maxwell's equations and applying the paraxial approximation we immediately get the time-domain wave equation where D = ε 0 E + P is the displacement field, and P is the induced polarization. In frequency domain, where E(r, ω) = ∞ −∞ dte iωt E(r,t) and similarly for P(ω), we may separate the linear and nonlinear polarization response as where χ (1) (ω) is the linear susceptibility tensor. The nonlinear part is usually expanded in a power series as The tensor products are calculated by summation over the indices ε −1 0 P NL,i = ∑ jk χ (2) i jk E j E k + ∑ jkl χ (3) i jkl E j E k E l + · · · . By assuming plane waves (no diffraction), we have where ε(ω) = 1 + χ (1) (ω) is the relative permittivity tensor, and we used ε 0 µ 0 = 1/c 2 . For SHG we define the electric field envelopes E j for the FW (frequency ω 1 ) and the SH (frequency ω 2 = 2ω 1 ), as well as envelopes for the nonlinear polarization response P NL, j E(z,t) = 1 2 u 1 E 1 (z,t)e ik 1 (ω 1 )z−iω 1 t + u 2 E 2 (z,t)e ik 2 (ω 2 )z−iω 2 t + c.c.
where the wave numbers are defined as k 2 j (ω) = ω 2 c 2 ε j (ω) with ε j (ω) = u j · χ (1) (ω), and where u is a unit polarization vector. In a lossless medium ε j is real, and the linear refractive index is given by n j (ω) = ε j (ω). The envelope wave equations are therefore The SVEA assumes | where we used the approximation k 2 j (ω) − k 2 j (ω j ) ≃ 2k j (ω j )[k j (ω) − k j (ω j )] [54, p. 34]. (A more accurate approach is to first go into the co-moving reference frame and then perform this expansion [55, p. 14], but eventually the SVEA will discard terms so the result will be the same.) At this time, we have also used that the assumption of a slow envelope compared to the fast oscillating carrier term e −iω j t means that the nonlinear polarization term on the right-hand side can be approximated as ω 2 2k j (ω j )c 2 ε −1 0 P NL, j (ω) ≃ ω 2 j 2k j (ω j )c 2 ε −1 0 P NL, j (ω). We can now perform the Taylor expansion of the wave number as k j (ω) = ∑ ∞ m=0 m! −1 k (m) j (ω m )(ω − ω j ), where the dispersion coefficients are defined as k (m) j (ω) = d m k j (ω)/dω m . Following this, we multiply both sides with e −iΩt where the frequency detuning is Ω = ω − ω j , and perform an inverse Fourier transform over Ω to time domain. This gives This is the SVEA equations in the stationary laboratory frame, describing a pulse envelope propagating with group velocity 1/k (1) j (ω j ). The time-domain dispersion operator iŝ and it can be truncated at any order to study effects higher-order dispersion. Besides a trivial transformation to a suitably chosen co-moving frame, what is left now is to evaluate the tensorial contributions to the nonlinear coefficients. Due to the vast variety of tensor combination possibilities for anisotropic nonlinearities it is here convenient to specify the problem very specifically: we therefore constrict ourselves to SHG in a uniaxial crystal, i.e. where light can be o-polarized with unit vector e o or e-polarized with unit vector e e [30, Eq. (B6)], and we also restrict ourselves to degenerate SHG, where the two FW photons are indistinguishable. Let us define the quadratic nonlinear polarization response as P (2) NL (z,t) = χ (2) : E(r,t)E(r,t), where we have assumed that the quadratic nonlinearity is instantaneous in time. The procedure is then to insert the field envelope in this nonlinear polarization response and then calculate the polarization response for a specific field by applying the dot product with the fields unit vector: ε −1 0 u j · P (2) NL (z,t). The only relevant terms are those that oscillate with the fast carrier frequency oscillations of the field e −iω j t , as this defines the polarization envelope P. Doing this we get for the FW and SH fields ε −1 0 P NL,1 (z,t) = χ (2) eff E * 1 (z,t)E 2 (z,t)e i∆kz , ε −1 0 P NL,2 (z,t) = χ where the phase mismatch is ∆k = k 2 (ω 2 ) − 2k 1 (ω 1 ). The effective susceptibility is calculated from the combination χ (2) eff = u j · χ (2) : u k u l that appears from isolating the e −iω j t term. The χ (2) is a rank 3 tensor with 27 different elements, but due to permutation symmetry in the the crystal axes indices there are only 18 nonzero elements [33]. Depending on the crystal symmetry class (point group) many of these are zero. It is therefore fruitful to consider a specific point group, and perform the calculations. We show in Sec. 3 the results for the point group 3m. We note that the nonlinear polarizations would have been a factor of 2 larger had we not used a factor 1 2 in front of our envelope definition.
Moving to the cubic nonlinear tensor component, it is done in much the same way, and we refer to [30, App. B] for details. We again consider an instantaneous response, and the cubic nonlinear polarization response that oscillates with e −iω j t turns out to be related to SPM and XPM; we could have had contributions from third-harmonic generation here if we had included such a harmonic field E 3 (z,t)e ik 3 (ω 3 )z−iω 3 t in the envelope definition. We then get [30, Eq. (B8)] ε −1 0 P NL, j (z,t) =

4 3χ
(3) eff (ω j ; ω j )|E j (z,t)| 2 + 6χ (3) eff (ω j ; ω k )|E k (z,t)| 2 E j (z,t) (47) The notation of the cubic susceptibilities χ eff (ω i ; ω j ) implicitly assumes phase-matched SPM and XPM interaction between ω i and ω j , as opposed to general four-wave mixing. We always have χ  1 2 , as in e.g. Boyd [28], it is important to note that the χ (3) eff values (and consequently also the n I 2,Kerr -values) in the two cases remain the same. The χ (3) is a rank 4 tensor with 81 different elements, but due to permutation symmetry in the the crystal axes indices there are only 30 nonzero elements [33]. The calculations of the χ We can now express the plane-wave SVEA equations for the electric field envelopes as eff (ω 2 ; ω 1 )E 2 |E 1 | 2 = 0 (49)