Second Harmonic Revisited : An Analytic Quantum Approach

We address the second-harmonic generation process in a quantum frame. Starting with a perturbative approach, we show that it is possible to achieve a number of analytic results, ranging from the up-conversion probability to the statistical properties of the generated light. In particular, the moments and the correlations of the photon-number distribution of the second-harmonic light generated by any initial state are retrieved. When possible, a comparison with the results achieved with the classical regime is successfully provided. The nonclassicality of some benchmark states is investigated by inspecting the corresponding autocorrelation function.


Introduction
Following the first observation of second-harmonic generation (SHG) by Franken et al. [1], the process has been the subject of a huge amount of theoretical and experimental investigations.The exact solution of the equation of motion for the classical field amplitudes, provided by Armstrong et al. [2], paved the way toward a formal analytic treatment of the problem.On the other hand, the quantum dynamics does not yield exact solutions, so that the quantum description of the process has been approached in a number of different ways, such as the analogy with the Dicke model [3,4] or the assumption of a factorization hypothesis on the moments of the photon number [5,6].However, the most common strategy is the perturbative approach [7], often followed by a numeric analysis [8][9][10].
Here, we review the simplest quantum description starting from a perturbative approach to the basic interaction Hamiltonian of the process, and we find a series of analytic results on the evolution of the pump field, the SHG probability, and the transformation of the statistics.In Section 2.1, we briefly report the well-known classical description of SHG focusing on the transformation of the statistics in the cases of coherent and chaotic light.In particular, we resume the results of a previous work by us, where, by applying a semiclassical approach, we have demonstrated [11] that the second-harmonic of a fundamental beam featuring multi-mode thermal statistics exhibits a more-than-thermal, namely superthermal, statistics.In Section 2.2, we approach the quantum model.Firstly, we retrieve some general results on the generation of second-harmonic (SH) light by inspecting the perturbative evolution of both the photon-number operator and the SH state.Then, we find that, in some simple cases, the SH conversion probability converges to periodic analytic functions.Finally, supported by the comparison with the classical results and by the quantum description, we move to the analysis of the output statistics by inspecting the Glauber autocorrelation function [12] of the SH field.In particular, we analytically retrieve it for the cases of coherent, chaotic, and squeezed light.

Classical Regime
Given an electric field E(x, t) propagating in a second-order nonlinear material, the SHG process can be described by solving Maxwell's equations in the absence of free charges and magnetic fields.The bound charges in the material are sources of a polarization vector P NL , which, in the particular case of SHG, is given by [2]: being 0 the permittivity of free space and χ (2) (ω SH = 2ω F ) the nonlinear susceptibility tensor related to this process.From now on, we fix the geometry and the symmetries of the nonlinear crystal so that χ (2) can be assumed to be a real scalar factor [13].In particular, we conform to the contracted notation and rename d eff ≡ χ (2) /2.Maxwell's equations provide a wave equation driven by this nonlinear polarization, i.e., [2,13]: where is the dimensionless permittivity of the material and c is the speed of light in free space.We assume that the propagation of the field is along the z axis and that the angular wave number is k j = n j ω j /c, where n j is the refractive index of the material corresponding to frequency ω j , being j = F, SH.
We express the fundamental field E F and the second-harmonic field E SH as [13]: and the nonlinear polarization as: NL (z, t) + P where [13] P It is common to assume the variation of the field amplitude small with respect to the field wavelength (slowly-varying-amplitude approximation), i.e., This assumption leads to the coupled-amplitude equations: where ∆k = 2k F − k SH is the phase mismatch.In the so-called undepleted-pump regime, we further assume that A SH (z) A F (z) ∀z, so that A F can be assumed to be constant with z, and the only equation to be solved in the system (8)-( 9) is the second one.
Hence, it is possible to retrieve the output intensity I SH = 2n SH 0 c|A SH | 2 as [13]: where L is the crystal length.
Given the link between the incoming light intensity and the second-harmonic one in Equation (10), it is possible to find a relation between the two corresponding light distributions [11,14,15].Upon writing Equation (10) as , where a is a coupling constant and f (I F ) is an invertible function (I F > 0), the second-harmonic light distribution G(I SH ) can be expressed as a function of the incoming one, P(I F ), which implies:

Coherent State
If the second-harmonic crystal is pumped with a coherent state, P(I F ) is a delta distribution: where I 0 is the intensity of the fundamental field at the crystal entrance.The second-harmonic distribution reads: By applying the properties of the delta distribution, Equation ( 14) takes the more familiar form: Therefore, the second-harmonic distribution in this case is the same as the input light distribution, which means that the detected-photon statistics of the incoming light and that of the second harmonic are both Poissonian.

Multithermal State
If the input light is a multithermal radiation with µ number of equally-populated modes, i.e., we find that the second-harmonic light displays superthermal fluctuations, according to: The first and second moments of the distribution in Equation ( 17) are: and the variance reads: Mandel's formula [16] provides the detected-photon distribution for the second-harmonic field [11], that is: where m and µ are the mean number of detected photons and the number of modes in the fundamental beam, Γ[j] is the Gamma function, and U[i, j, k] is the Tricomi confluent hypergeometric function.Upon introducing the quantum detection efficiency ζ, the mean value and the variance of this distribution are: We recall that the convolution of a large number of thermal distributions converges to a Poissonian distribution in the discrete case (i.e., the detected-photon statistics) and to a delta distribution in the continuous case (i.e., the light intensity statistics).It is worth noting that this is true also for a superthermal distribution, as we can see from Equations ( 20) and (23) since: In particular, the convergence of the superthermal distribution to a Poissonian is four-times faster than that of the multithermal distribution.In fact, for µ 1:

Quantum Regime
The quantum dynamics of second-harmonic generation can be described through the Hamiltonian: where âj and â † j are the annihilation and creation boson operators for the input (j = F) and second-harmonic field (j = SH) with the usual commutation rules âj , â † k = δ j,k .No exact solution can be found for the output state since it is not possible to identify a finite-dimensional Lie algebra for this Hamiltonian [17].As mentioned in the Introduction, here we exploit a perturbative approach.A comparison with the classical regime will help to test the consistency of our results.
The Heisenberg equation for the number operator nSH ≡ â † SH âSH , with the Hamiltonian Ĥ in Equation (25), has no exact solution, as mentioned above.Nevertheless, it allows one to retrieve the coefficients of the Taylor expansion for nSH , being the derivatives of all orders provided by the iteration of Equation ( 26) and n0 ≡ nSH (0).
Equivalently, one can study the evolution of the states.The most general input state for an SHG process can be expanded over the basis of the number states as: where c k are suitable coefficients such that Given the evolution operator Û = exp i h Ĥt , we get the output state from: by suitably expanding Û.In particular, we focus on the evolution of Fock states and then generalize to any input state by expanding it like in Equation (28), i.e., Thus, we just need to compute |ψ(γt) = Û † |ψ(0) , where |ψ(0) ≡ |N, 0 .For instance, up to the second order in γ, the output state reads: Both the Heisenberg (Equation ( 27)) and the Schrödinger (Equation ( 29)) pictures will be exploited in the following.In particular, in Section 2.2.1, we need the expansion of the number operator to retrieve the first moment of the SH photon-number distribution.On the contrary, the conversion probabilities in Section 2.2.2 are more easily computed by evolving the input state.

Photon-Number Distribution
Kozierowski and Tanaś in [7] provided an expansion of the first moment of the output statistics â † SH âSH up to the fourth order as a function of the Glauber autocorrelation function of the input field.The Glauber autocorrelation function for the modes âF and â † F is defined as follows [12,16]: As shown in Appendix A, it is possible to demonstrate that every order of the perturbative expansion of the output mean photon-number can be expressed as a linear combination of the Glauber autocorrelation functions in Equation (32).
The theorem allows us to express the SH moments as convergent series over the order of the autocorrelation function of the fundamental field.In particular, we found out that up to the sixth order in γ, the output mean value reads: F N 2 + 2g (3) where N ≡ â † F âF .For the sake of simplicity, we assumed γ ∈ R and replaced γt with γ.
Moreover, for the Fock state, the theorem provides a general expression for the asymptotic value of nSH as a function of N (see Appendix B), namely: Note that, upon defining the conversion efficiency as: we find from Equation (34) that: meaning that the conversion efficiency related to every single interaction is proportional to the input intensity.This fact will be relevant in the analysis of the evolution of the input light statistics.Furthermore, note that, if γ2 N < 1, Equation ( 34) is a geometric series, and nSH ∼ γ2 N 2 /(1 − γ2 N).
The physical meaning of the coefficients in the perturbative terms of nSH can be pointed out by writing the autocorrelations in Equation (33) for the Fock state as a function of N: for any given K. Up to the sixth order in γ, we have: Note that the three terms in Equation (37) represent the three different processes concurring in SHG.The first one is related to the creation of k SH photons from N. This contribution comes from all and only the terms that in the expansion of U in Equation ( 29) are of the kind ( â † SH ) k .They correspond to pure up-conversion events, i.e., the creation of 2k SH photons from N at each order k.
Then, for k ≥ 3, we have the α k terms, corresponding to the alternating creation and annihilation of SH photons.For the sake of simplicity, we name these as annihilation processes since, at variance with the former term, they consist of back-conversion events.For instance, for k = 3, we have three different processes leading to a non-null value of nSH , which are the following: the creation process ( â † SH ) 3 (creation of three SH photons from N after three interactions) and the annihilation processes â † SH âSH â † SH and âSH â † SH â † SH (creation of one SH photon from N after three interactions).The last two provide the first and the second term in α 3 (Equation (38)), respectively, and differ for the ordering of the operators only.In particular, both these processes lead to the same output state (|N − 2, 1 ), but, in the first one, we have that an SH photon is created, annihilated, and then created again (each of these events happening in N!/(N − 2)! possible ways), while in the second one, two SH photons are created (which is likely to happen in N!/(N − 4)! equivalent configurations), and then, one of the two is down-converted in two photons of the fundamental field ((N − 2)!/(N − 4)! possible configurations).
Finally, if k ≥ 2, the annihilation processes can interfere either with one another or with creation events, giving rise to the third kind of contribution (the β k terms in Equation ( 37)).This superposition happens whenever different processes lead to the same output state.For instance, for k = 2, we find, together with the two-photon-creation process ( â † SH ) 2 (creation of two SH photons from N after two interactions), the superposition β 2 of the one-photon-creation process â † SH (creation of one SH photon from N after one interaction) with the annihilation processes â † The role of these different contributions will be outlined in the following, in connection with the conversion probability and the evolution of the input light statistics.

Conversion Probability
Achieving an analytic expression for the probability of converting 2k photons of the fundamental field into l ≤ k SH photons is generally a hard task because the number of perturbative terms in the expansion of the SH state (see Equation ( 29)) is generally not enough to recognize the convergence to a known analytic function.Moreover, there is no guarantee at all that such a function exists.However, there are some simple cases where a large number of perturbative orders (up to 30) can be retrieved and found to be the terms of the perturbative expansion of known analytical functions.In particular, this is the case if the input photon-number is small and l = 1.We found out that the probability p 1 (N, γ) of generating a single SH photon from 2 ≤ N < 6 input photons converges, up to the order 30, to a periodic function, i.e., where: The interplay among multiple creation and annihilation processes and their superposition is periodic up to six input photons, as shown in Figure 1.For N ≥ 6, this is no longer true, but still, p 1 (N, γ) converges to an analytic function.The reason for the aperiodic behavior of p 1 (N, γ) for N ≥ 6 is due to the contribution of a larger number of processes at every order, resulting in a squared sum of different periodic functions.The ratio between the periods of these functions turns out to be incommensurable, so that their sum is not a periodic function.

Statistics
Finally, we focus on the transformation of the statistics of the input field via SHG.We choose a particular input statistics by suitably setting the coefficients c n of the input density matrix (28).We investigate and compare the evolution of Fock, coherent, multithermal, and squeezed states.Moreover, we point out the role of the creation and annihilation processes together with their superpositions, outlined in Equation (38) by comparing the autocorrelation function for SHG (blue solid lines in the plots) with the autocorrelation function for the process ∑ K k=1 γ2k N! k!(N−2k)!k (red dashed lines in the plots), only corresponding to the SH creation events.The latter is simply obtained by neglecting the α k and β k terms in Equation (37).Since the conversion efficiency is proportional to the input photon-number (see Equation ( 36)), we expect the annihilation processes and the superpositions to reduce the quantum correlations in the output state.The moments of the output statistics are straightforward from: In analogy with Equation (33), the variance can be expressed as an expansion over the autocorrelation functions (32) of the input field, which here we rename g F : F N 2 + 2g (3) F N 2 + 2g (3) It is possible to show that the autocorrelation function defined in Equation (32) can be written as a function of the variance and of the mean value of the number operator, so that we can compute the autocorrelation function of the SH field through:

Fock State
From Equation (41), we immediately get: which is consistent with the theoretical expectations since it is a fact (see Appendix C) that, up to the first order, g The second-and fourth-order autocorrelation functions for a Fock state respectively read: so that, applying Equation ( 43), the first-order expansion for the expected second-harmonic autocorrelation function must be:

Coherent State
The case of a coherent input state can be approached by studying the evolution of Equation (28) upon requiring: where |α| is the amplitude of the generic coherent state |α = e −|α| 2 ∑ ∞ k=0 |α| 2k /k! |k .
Up to k = 2, the autocorrelation function reads: and shows that the second-harmonic output state displays a sub-Poissonian statistics at the first perturbative order.Thus, the SH field from a coherent pump field is endowed with quantum correlations, which is due to the anti-bunching resulting from an up-conversion with efficiency proportional to N j (Equation (36)).It is a well-known result [8,9], which has been also experimentally tested [18,19].We here point out that, as expected, the annihilation processes and the interference between different processes contribute to reducing quantum correlations and, consequently, the nonclassicality of the output state, as shown in Figure 2.

Multithermal State
The probability distribution of a multithermal state reads [16]:

Discussion
We investigated the second-harmonic generation process starting from the usual perturbative approach, but aiming to achieve some analytic results.In particular, we focused on the statistics of the output field, and helped by the comparison with the well-known classical regime, we retrieved the autocorrelation functions of the SH field for different input states, which are, as far as we know, of experimental interest.In fact, among the nonlinear parametric processes, the SHG is widely explored also in view of its multiple applications, which range from the production of laser systems with multiple wavelengths to the matching of the operation spectral range of photodetectors [21].Moreover, nonclassical states of light are an essential resource for quantum technologies.Thus, investigating the way in which quantum optical states can be obtained through the SHG process represents a fundamental step towards their tailoring and exploitation.For instance, only a quantum approach predicts the generation of a sub-Poissonian state starting from a coherent one, a fact that is also experimentally observable.
In addition, we tried to reach a deeper understanding of the SHG process by approaching the quantum problem with different strategies and looked for the physical meaning of the perturbative terms in the expansion of the moments.In particular, we found that, at least in some cases, the conversion probability can be written in a closed form.Albeit that these latter results are not amenable for a direct experimental test, we think that they could be useful for further analytic predictions.On the one hand, an analytic formula up to a generic order k for the terms α k and β k in Equation (37), specified in Equation (38) for the first three orders, would allow us to retrieve a general expression for the moments of the SH statistics, and on the other hand, further expansions of the output SH state will lead to computing the conversion probability in more complicated cases.These two issues are both non-trivial since they require a larger number of perturbative orders, which is a computationally hard task.every even derivative, but lowered by one for every odd derivative, so that k ≤ j 2 + 1.Hence, we can conclude that n is a linear combination of ( â † F ) k âk F = γk g (k) F N k .

Appendix B
It is well known that the autocorrelation function of order k for a Fock state |N reads: As mentioned in Section 2.2, the transformation of the number state is the starting point for the analysis of the evolution of the other light distributions under investigation, from which the search of the implications of the theorem in Appendix A for an input number state follows.
Corollary A1.If the fundamental input mode is a Fock state |ψ(0) = |N, 0 , then: Proof. ( â † ) k âk = N!/(N − k)! ∼ N k .The highest k at the order j is j 2 + 1, but, since the odd j terms do not contribute, we only sum over the even j orders and replace j with 2j.

Appendix C
Up to the first order, the autocorrelation function can be approximated to Equation (43).The result is obtained through a first-order assumption on the modes, i.e., a SH = γa F a F .Then, from Equation (41), we find:

Figure 1 .
Figure 1.Evolution of the probability of generating one SH photon from N pump photons.Solid red line: Tr[ ρ( γ)1 F ⊗ nk ] where ρ is obtained by expanding Û in Equation (29) up to k = 30.Dashed line: theoretical prediction from Equation (39).The purple solid line in the bottom left panel is the probability of generating two photons from four.It is of the form A sin 4 (ω γ).

Figure 2 .
Figure 2. Coherent input state.Solid blue line: SH autocorrelation function, from Equation (46).Dashed line: SH autocorrelation function in the absence of annihilation processes and superpositions.