1 Introduction

The application of Kramers–Kronig (K–K) dispersion relations in data inversion of linear optical constants has been widely used since the 1920s [1], which allows the calculation of the linear refractive index from the measured absorption spectra and vice versa. However, it took several decades until the rigorous proof of the equivalence between causality and the existence of dispersion relations was given by Toll [2]. The principle of causality guarantees the holomorphicity (or analyticity) of the linear optical constant i.e. the poles (= singular points of the function) are located at the lower half of the complex frequency plane allowing the derivation of the conventional K–K relations [3]. A review that summarized well the connection between causality and dispersion relations was published recently [4] that extended the analysis also to metamaterials.

Utilization of dispersion relations requires the knowledge of the measured part over the complete spectral range since the dispersion integrals are calculated from zero to infinity. It has been shown that extrapolations of data beyond the measured spectral range can be a source of serious errors and dispersion relation should be used with care [5, 6]. Therefore, it is desirable to search for dispersion relations with fast convergence. On one hand, the convergence of K–K relations can be improved if one has a priori knowledge of the unknown part at some distinct frequencies. These so called anchor points allow the derivation of singly (SSKK) [7] or multiply (MSKK) [8] subtractive K–K relations, which have a faster convergence than the conventional K–K relations. On the other hand, by considering the powers of the optical constants [9] one observes even faster convergence than with the subtractive relations.

The invention of laser resulted in the fast development of nonlinear optics [10] and the first extensions of dispersion relations for nonlinear optical susceptibilities took place already in the 1960s [11,12,13]. However, a detailed investigation of the asymptotic behaviour and the holomorphicity of nonlinear susceptibilities was not started until the 1990s by Bassani et al. [14, 15]. Peiponen [16] was the first to study nonlinear sum rules, which are general constraints that either model generated or measured data must obey, and they give information of the optical response. Unfortunately, when dealing with nonlinear susceptibilities, the susceptibility may be a meromorphic function i.e. the function may have poles in both halves of the complex frequency plane. Hence, the conventional K–K relations are invalid and the effect of the poles in the upper half plane must be taken into account by the theorem of residues [17]. Such generalized K–K type relations and sum rules have been derived for pump and probe experiments [18] and for nonlinear refractive index [19]. For nonlinear susceptibilities the convergence of the K–K type relations plays a crucial role since the asymptotic fall-off of the nonlinear susceptibility increases with the order of the nonlinear process [20]. On contrary, the signal of a nonlinear process decreases with the order of the nonlinear process [21]. Therefore, the extension of dispersion theory for the moments and the powers of meromorphic optical susceptibilities exploiting strong convergence of dispersion integrals is important in testing of theoretical models describing the nonlinear response.

Nowadays the spectral properties of the nonlinear refractive index and two-photon absorption are important, for example, in the development of new optoelectronic devices based on nanostructures such as quantum wire lasers and all-optical switching devices [22], and in modern drug research based on the use of nanoparticles for bioaffinity assays [23]. The latter technique has already led to a portable commercial device that can be used for a relatively fast detection of hospital bacteria such as methicillin-resistant Staphylococcus aureus (MRSA). Another field of application of the present dispersion relations is related to THz generation via four-wave mixing in plasmas, where the third-order nonlinear susceptibility plays an important role [24].

Generalized Kramers–Kronig relations for meromorphic nonlinear susceptibility have been presented previously in our book [17]. Here the theory is extended by considering generalized dispersion relations and sum rules for the moments and powers of degenerate four wave mixing susceptibility, which is related to the nonlinear refractive index and two-photon absorption. Such relations and sum rules have not been presented until now.

Recently the degenerate four wave mixing susceptibility has been studied in various nanostructures ranging from nanodisks [25] and graphene [26] to doped thin films [27] in which the nonlinear signal can be significantly enhanced due to high local electric fields. The generalized dispersion relation formalism developed in this paper can be applied for a more accurate analysis of such enhanced nonlinear response providing insight and validation of the degenerate nonlinear four wave mixing susceptibility.

In Sect. 2 a general response function of the degenerate third-order susceptibility is derived that is a meromorphic function. Furthermore, a classical anharmonic oscillator model is applied to obtain the asymptotic fall-off and symmetry properties of the moments of the degenerate third-order susceptibility, and the model is further used in the calculation of the analytic expressions for the poles in the upper half plane. The classical K–K relations are briefly derived and presented in Sect. 3. In Sect. 4 theorem of residues is applied for the poles in the upper half plane followed by generalized K–K type relations and sum rules. Is is noteworthy that the results presented in Sect. 4 are general and applicable to any meromorphic optical constant. Finally, Sect. 5 summarizes the findings.

2 General properties of the degenerate four wave mixing susceptibility

The electrical polarization P(t) generated in a media by the presence of an electric field E(t) to an arbitrary-order n in perturbation theory is controlled by the response function \(G^{(n)}\). The polarization of a media up to third-order is given by

$$\begin{aligned} P_i(t)= & {} \int _0^{\infty }G_{ij}^{(1)}(t_1)E_j(t-t_1){\mathrm{d}}t_1\nonumber \\&+\int \int _0^{\infty }G_{ijk}^{(2)}(t_1,t_2)E_j(t-t_1)E_k(t-t_2){\mathrm{d}}t_2{\mathrm{d}}t_1\nonumber \\&+\int \int \int _0^{\infty }G_{ijkl}^{(3)}(t_1,t_2,t_3)E_j(t \nonumber \\&-t_1)E_k(t-t_2)E_l(t-t_3){\mathrm{d}}t_3{\mathrm{d}}t_2{\mathrm{d}}t_1, \end{aligned}$$
(1)

where indices ijkl denote the Cartesian components and G’s are real valued response functions. The response of a medium can only depend on preceding times and therefore, the response function must be equal to zero for all times \(t_{1,2,3}<0\) (the principle of causality). In this paper the third-order response in degenerate four wave mixing is studied. A Fourier-transform of the third-order response function gives the general third-order susceptibility as follows:

$$\begin{aligned} \chi _{ijkl}^{(3)}(\omega _1,\omega _2,\omega _3)= & {} \int \int \int _0^{\infty }G_{ijkl}^{(3)}(t_1,t_2,t_3) \exp [{\mathrm{i}}(\omega _1t_1+\omega _2t_2 \nonumber \\&+\omega _3t_3)]{\mathrm{d}}t_3{\mathrm{d}}t_2{\mathrm{d}}t_1. \end{aligned}$$
(2)

In the 3-dimensional frequency space one can define, to obtain dispersion relations and sum rules, a general line of integration [14]. Here, the degenerate four wave mixing is considered i.e. the interaction of three coherent light beams with the same frequency. Therefore, the third-order response as a function of \({\varOmega }\) is of the form

$$\begin{aligned} \chi _{ijkl}^{(3)}({\varOmega },-{\varOmega },{\varOmega })= & {} \int \int \int _0^{\infty }G_{ijkl}^{(3)}(t_1,t_2,t_3) \exp [{\mathrm{i}}(t_1-t_2 \nonumber \\&+t_3){\varOmega }]{\mathrm{d}}t_3{\mathrm{d}}t_2{\mathrm{d}}t_1 \end{aligned}$$
(3)

Unfortunately, the analytic continuation of the above susceptibility in the complex \({\varOmega }\)-plane is a meromorphic function [14], and the conventional K–K relations are invalid.

As an example a classical anharmonic oscillator is used that gives the correct asymptotic fall-off and symmetry properties of the degenerate third-order susceptibility. The general nature of dispersion relations should be emphasized here as they depend only on the complex nature and the asymptotic fall-off of the investigated function. Moreover, it should be noted that the results obtained in the next section are valid not only in the frame of the anharmonic oscillator model but for all meromorphic functions, which have the proper asymptotic fall-off and symmetry properties. The anharmonic oscillator model for degenerate four wave mixing susceptibility is given by [28]

$$\begin{aligned} \chi ^{(3)}({\varOmega },-{\varOmega },{\varOmega })=B[D({\varOmega })]^3D(-{\varOmega })=B|D({\varOmega })|^2[D({\varOmega })]^2, \end{aligned}$$
(4)

where B is a constant, which gives information about electron number density and strength of the anharmonic terms in equation of motion of the electron, \({\varOmega }\) is the complex frequency and \(D({\varOmega })\) is given by

$$\begin{aligned} D({\varOmega })=(\omega _0^2-{\varOmega }^2-{\mathrm{i}}{\varGamma }{\varOmega })^{-1}. \end{aligned}$$
(5)

Here \(\omega _0\) and \({\varGamma }\) are the resonance frequency and the half width of the resonance peak, respectively. For the sake of simplicity the Cartesian indices have been dropped from the following notation. The third-order poles of \([D({\varOmega })]^3\) are located in the lower half plane. The meromorphicity of \(\chi ^{(3)}({\varOmega },-{\varOmega },{\varOmega })\) is caused by a term \(D(-{\varOmega })\) in Eq. (4), which has two first-order poles in the upper half plane. Therefore, the degenerate third-order susceptibility is a meromorphic function both in the lower and the upper half planes. Furthermore, the real part of the degenerate third-order susceptibility is connected to the nonlinear refractive index \(n_{\mathrm{{NL}}}({\varOmega })\) whereas the imaginary part of the susceptibility is related to the two-photon absorption \(\kappa _{\mathrm{{NL}}}({\varOmega })\) [21]. The generalized K–K relations obtained in Sect. 4 for the degenerate third-order susceptibility can thus be used to evaluate both the nonlinear refractive index and the two-photon absorption coefficient.

In this paper the results given by Bassani and Lucarini [18] are generalized by considering the moments and the powers of a third-order degenerate susceptibility, which have very rapid convergence. Hence, the errors in data inversion caused by finite spectral range are assumed to diminish. This assertion has been demonstrated numerically using the measured data of holomorphic third-harmonic generation susceptibility of polymers polytiophene and polysilane [29]. It was observed that the upper limits for the moments and the powers of the third-order harmonic generation are valid and the higher moments are divergent as predicted in Ref. [30].

The asymptotic behaviour of the third-order degenerate susceptibility can be derived either quantum mechanically or from the anharmonic oscillator model, which both yield

$$\begin{aligned} \chi ^{(3)}({\varOmega },-{\varOmega },{\varOmega })=\psi {\varOmega }^{-8}+O({\varOmega }^{-8}), \end{aligned}$$
(6)

where \(\psi\) is a high frequency constant of the susceptibility and \(O({\varOmega }^{-8})\) denotes terms that converge strictly faster than \({\varOmega }^{-8}\). The requirement of the real valued response functions leads to a symmetry property

$$\begin{aligned} \chi ^{(3)}(-{\varOmega }^{*},{\varOmega }^{*},-{\varOmega }^{*})=[\chi ^{(3)}({\varOmega },-{\varOmega },{\varOmega })]^{*} \end{aligned}$$
(7)

with \((^{*})\) denoting complex conjugation. Hence, the real part of the nonlinear susceptibility is an even function and the imaginary part an odd function.

3 K–K relations and sum rules for holomorphic linear optical constants

First a response function of a linear susceptibility \(\chi ^{(1)}(\omega )\) is considered that is given by

$$\begin{aligned} \chi ^{(1)}({\varOmega })=\int _0^{\infty }G^{(1)}(t) \exp ({\mathrm{i}}{\varOmega }{t}){\mathrm{d}}t. \end{aligned}$$
(8)

A simple solution for the response function is given by the classical Lorentz oscillator model with a plasma frequency of \(\omega _{\mathrm{{p}}}^2\) as follows [31]:

$$\begin{aligned} \chi ^{(1)}({\varOmega })=-\omega _{\mathrm{{p}}}^2D({\varOmega }) \end{aligned}$$
(9)

with \(D({\varOmega })\) given by Eq. (5). The asymptotic behaviour of the linear complex susceptibility of the Lorentz model is proportional to \(|{\varOmega }|^{-2}\) at high angular frequencies. Moreover, \(D({\varOmega })\) has two complex poles that are located in the lower half plane and are explicitly given by

$$\begin{aligned} {\varOmega }_{1,2}=\pm \sqrt{\omega _0^2-{\varGamma }^2/4}-{\mathrm{i}}{\varGamma }/2. \end{aligned}$$
(10)

Therefore, the Lorentzian susceptibility \(\chi ^{(1)}(\omega )={{\mathrm{Re}}}\{\chi ^{(1)}(\omega )\}+{\mathrm{i}}{{\mathrm{Im}}}\{\chi ^{(1)}(\omega )\}\) is a holomorphic function in the upper half plane and the K–K relations can be written as follows:

$$\begin{aligned} {{\mathrm{Re}}}\{\chi ^{(1)}(\omega ')\}=\frac{2}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{\omega {{\mathrm{Im}}}\{\chi ^{(1)}(\omega )\}}{\omega ^2-\omega '^2}{\mathrm{d}}\omega , \end{aligned}$$
(11)
$$\begin{aligned} {{\mathrm{Im}}}\{\chi ^{(1)}(\omega ')\}=-\frac{2\omega '}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{{{\mathrm{Re}}}\{\chi ^{(1)}(\omega )\}}{\omega ^2-\omega '^2}{\mathrm{d}}\omega . \end{aligned}$$
(12)

Kramers [32] used the definition of complex refractive index \(N(\omega )=n(\omega )+{\mathrm{i}}k(\omega )\), where \(n(\omega )\) is the frequency dependent real refractive index of a medium and \(k(\omega )\) the extinction coefficient and showed that the complex refractive index obeys the symmetry relation \(N(-\omega )=[N(\omega )]^{*}\). The asymptotic fall-off of the complex refractive index is proportional to \(|{\varOmega }|^{-2}\) [33], similar to the linear dielectric function. Hence, the K–K relations hold true for the complex refractive index and are given by

$$\begin{aligned} n(\omega ')-1=\frac{2}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{\omega k(\omega )}{\omega ^2-\omega '^2}{\mathrm{d}}\omega , \end{aligned}$$
(13)
$$\begin{aligned} k(\omega ')=-\frac{2\omega '}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{n(\omega )-1}{\omega ^2-\omega '^2}{\mathrm{d}}\omega . \end{aligned}$$
(14)

Equations (13) and (14) are the classical Kramers–Kronig relations. Classical K–K relations take into account only the frequency dependent optical constants of the medium. With the aid of the K–K relations, one can calculate the refractive index of a medium from measured absorbance spectrum and vice versa.

The importance of sum rules in linear optical spectroscopy should not be underestimated. Perhaps the best known sum rule is the f-sum rule that can be derived from a K–K relation between the real and imaginary parts of linear complex susceptibility, \(\chi ^{(1)}(\omega )={{\mathrm{Re}}}\{\chi ^{(1)}(\omega )\}+{\mathrm{i}}{{\mathrm{Im}}}\{\chi ^{(1)}(\omega )\}\). The f-sum rule is as follows:

$$\begin{aligned} \int _0^{\infty }\omega {{\mathrm{Im}}}\{[\chi ^{(1)}(\omega )]\}{\mathrm{d}}\omega =\frac{\pi }{2}\omega _{\mathrm{{p}}}^2, \end{aligned}$$
(15)

where \(\omega _{\mathrm{{p}}}\) is the plasma frequency and provides important information about electron number density of the medium. One may consider that integration over semi-infinite frequency axis is an issue in practical use of the f-sum rule. Due to the access of optical constants in rather wide spectral range this issue has become less severe as clearly pointed in a recent study of Franta et al. [34]. Another sum rule that is obtained from a K–K relation is a static limit sum rule between refractive index \(n(\omega )\) and extinction coefficient \(k(\omega )\). This relation is given by the expression

$$\begin{aligned} n(0)-1=\frac{2}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{k(\omega )}{\omega }{\mathrm{d}}\omega , \end{aligned}$$
(16)

A recent study of optical properties of water vapor in terahertz spectroscopy by Grischkowsky et al. [35] emphasized the importance of \(n(0)-1\) in correct description of the spectral properties of water vapor.

4 Generalized meromorphic K–K type relations and sum rules

Fig. 1
figure 1

Contour for derivation of Kramers–Kronig type dispersion relations and sum rules (\(\times\)=pole)

Full size image

Dispersion theory for meromorphic optical constants is based on the utilization of complex contour integration along the closed integration path C presented in Fig. 1. Three poles are located in the upper half plane: two of them are, due to the meromorphic nature of the degenerate third-order susceptibility, located symmetrically with respect to the imaginary axis, and one pole is on the real axis due to the definition of the function under investigation i.e.

$$\begin{aligned} f({\varOmega })=\frac{{\varOmega }^m[\chi ^{(3)}({\varOmega },-{\varOmega },{\varOmega })]^k}{{\varOmega }-\omega '}, \end{aligned}$$
(17)

where \(\omega '\) is a real valued frequency and m and k are positive integers. The asymptotic fall-off is obtained from Eq. (6) and symmetry property is given by Eq. (7). To obtain dispersion relation, a complex integration of the function \(f({\varOmega })\) is carried out along the path presented in Fig. 1. The contour integration along the closed circle C can be split into parts after the limiting process \(|{\varOmega }|\rightarrow \infty , \rho \rightarrow 0\) as follows:

$$\begin{aligned} \oint _{\mathrm{{C}}}f({\varOmega })\,{\mathrm{d}}{\varOmega }= & {} {\mathrm{P}}\int _{-\infty }^{\infty }f(\omega )\,{\mathrm{d}}\omega +\lim _{\rho \rightarrow 0}\int _{\gamma }f({\varOmega })\,{\mathrm{d}}{\varOmega } \nonumber \\&+\lim _{|{\varOmega }|\rightarrow \infty }\int _{\mathrm{{A}}}f({\varOmega })\,{\mathrm{d}}{\varOmega }, \end{aligned}$$
(18)

where P indicates Cauchy’s principal value integration and \(\gamma\) denotes the clockwise integration path along the small semi-circle whereas A is the counterclockwise integration path along the large semi-circle. Since the investigated function \(f({\varOmega })\) is meromorphic in the complete complex plane, the poles in the upper half plane need to be taken into account. According to Cauchy’s integral theorem and the theorem of residues [17], the complex integration around the closed circle C is followed by

$$\begin{aligned} \oint _{\mathrm{{C}}}f({\varOmega })\,{\mathrm{d}}{\varOmega }=2\pi {\mathrm{i}}\sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\,[f({\varOmega })]. \end{aligned}$$
(19)

Here, the summation of the residues is over the poles in the upper half plane. Similarly for the first-order pole at \({\varOmega }=\omega '\), which is located exactly at the integration path, one obtains

$$\begin{aligned} \lim _{\rho \rightarrow 0}\int _{\gamma }f({\varOmega })\,{\mathrm{d}}{\varOmega }=-{\mathrm{i}}\pi {\mathrm{Res}}_{{\varOmega }=\omega '}f({\varOmega })=\omega '^m[\chi ^{(3)}(\omega ',-\omega ',\omega ')]^k. \end{aligned}$$
(20)

The upper limits for the convergent moments m are obtained by the investigation of counterclockwise integration along the semi-circle A in the high frequency limit (for a detailed derivation see Appendix B in Ref. [30]).

Assuming that m is an even integer allows derivation of K–K type dispersion relations. For the sake of simplicity, a simple notation \(\chi ^{(3)}(\omega )\) will be used instead of the full \(\chi ^{(3)}(\omega ,-\omega ,\omega )\) in the following derivations. The symmetry property of Eq. (7) is followed by

$$\begin{aligned}&{\mathrm{P}}\int _{-\infty }^{\infty }\frac{\omega ^m[\chi ^{(3)}(\omega )]^k}{\omega -\omega '}\,{\mathrm{d}}\omega \nonumber \\&\quad =2{\mathrm{P}}\int _{0}^{\infty }\frac{\omega '\omega ^m{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega )]^k\}+{\mathrm{i}}\omega ^{m+1}{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega )]^k\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega . \end{aligned}$$
(21)

By splitting Eqs. (19), (20), and (21) into the real and imaginary parts, the following K–K type relations are obtained for even integers m with \(m\le 8k-2\) (this asymptotic limit is obtained from the convergence of integrands by a similar derivation given in Ref. [30])

$$\begin{aligned}&\omega '^m{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega ')]^k\}=\frac{2}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{\omega ^{m+1}{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega )]^k\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega \nonumber \\&\quad -2{{\mathrm{Re}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\left\{ \frac{{\varOmega }^m[\chi ^{(3)}({\varOmega })]^k}{{\varOmega }-\omega '}\right\} \right\} , \end{aligned}$$
(22)
$$\begin{aligned}&\omega '^m{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega ')]^k\}=-\frac{2\omega '}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{\omega ^m{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega )]^k\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega \nonumber \\&\quad -2{{\mathrm{Im}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\left\{ \frac{{\varOmega }^m[\chi ^{(3)}({\varOmega })]^k}{{\varOmega }-\omega '}\right\} \right\} . \end{aligned}$$
(23)

For the sake of simplicity, a simple notation for the degenerate third-order susceptibilty analogous to Eqs. (13) and (14) is obtained by setting \(k=m=0\) as follows:

$$\begin{aligned} {{\mathrm{Re}}}\{[\chi ^{(3)}(\omega ')]\}= & {} \frac{2}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{\omega {{\mathrm{Im}}}\{[\chi ^{(3)}(\omega )]\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega \nonumber \\&-2{{\mathrm{Re}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\left\{ \frac{[\chi ^{(3)}({\varOmega })]}{{\varOmega }-\omega '}\right\} \right\} , \end{aligned}$$
(24)
$$\begin{aligned} {{\mathrm{Im}}}\{[\chi ^{(3)}(\omega ')]\}= & {} -\frac{2\omega '}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega )]\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega \nonumber \\&-2{{\mathrm{Im}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\left\{ \frac{[\chi ^{(3)}({\varOmega })]}{{\varOmega }-\omega '}\right\} \right\} . \end{aligned}$$
(25)

In a similar way for an odd integer m the following relation is obtained

$$\begin{aligned}&{\mathrm{P}}\int _{-\infty }^{\infty }\frac{\omega ^m[\chi ^{(3)}(\omega )]^k}{\omega -\omega '}\,{\mathrm{d}}\omega \nonumber \\&\quad =2{\mathrm{P}}\int _{0}^{\infty }\frac{\omega ^{m+1}{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega )]^k\}+{\mathrm{i}}\omega '\omega ^m{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega )]^k\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega , \end{aligned}$$
(26)

which is followed by K–K type relations for an odd \(m\le 8k-1\) as follows:

$$\begin{aligned}&\omega '^m{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega ')]^k\}=\frac{2\omega '}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{\omega ^m{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega )]^k\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega \nonumber \\&\quad -2{{\mathrm{Re}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\left\{ \frac{{\varOmega }^m[\chi ^{(3)}({\varOmega })]^k}{{\varOmega }-\omega '}\right\} \right\} , \end{aligned}$$
(27)
$$\begin{aligned}&\omega '^m{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega ')]^k\}=-\frac{2}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{\omega ^{m+1}{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega )]^k\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega \nonumber \\&\quad -2{{\mathrm{Im}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\left\{ \frac{{\varOmega }^m[\chi ^{(3)}({\varOmega })]^k}{{\varOmega }-\omega '}\right\} \right\} . \end{aligned}$$
(28)

For the value \(m=8k\), which is the highest convergent moment, the complex integration along the semi-circle A does not converge to zero but yields information of the high frequency value of the degenerate susceptibility. For higher values m the integral diverges. Therefore, the K–K type relation (22) is transformed as follows:

$$\begin{aligned}&\omega '^m{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega ')]^k\}=\frac{2}{\pi }{\mathrm{P}}\int _0^{\infty }\frac{\omega ^{m+1}{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega )]^k\}}{\omega ^2-\omega '^2}\,{\mathrm{d}}\omega \nonumber \\&\quad -2{{\mathrm{Re}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}\frac{{\varOmega }^m[\chi ^{(3)}({\varOmega })]^k}{{\varOmega }-\omega '}\right\} +\psi ^k, \end{aligned}$$
(29)

whereas the K–K type relation (23) remains in the same form. The expressions (22), (23), (27) and (28) are generalizations of Kramers–Kronig relations in the sense that K–K relations are “corrected” by residue terms, and hold for powers and moments of the degenerate nonlinear susceptibility.

The poles of the degenerate third-order susceptibility are of the first-order as shown in Sect. 2. This was taken into account by Bassani and Lucarini [18] as they considered explicitly the poles in the upper half plane for pump and probe nonlinear processes in the frame of the anharmonic oscillator model, which provides exact results as concerns the residue terms. In this paper the previous results [18] are generalized and obtain better convergence with generalized K–K relations and sum rules by considering the powers of the degenerate nonlinear susceptibility. Unfortunately, for the poles of the order \(\ge 2\) the derivatives of the complex function appear. Therefore, since the kth power of the degenerate susceptibility has two kth-order poles in the upper half plane, denoted by \({\varOmega }_1\) and \({\varOmega }_2\), the residue terms result with poles higher order than 1 as follows:

$$\begin{aligned}&{\mathrm{Res}}_{{\varOmega }={\varOmega }_{1,2}}\frac{[\chi ^{(3)}({\varOmega })]^k}{{\varOmega }-\omega '} \nonumber \\&\quad =\lim _{{\varOmega }\rightarrow {\varOmega }_{1,2}}\frac{1}{(k-1)!}\frac{{\mathrm{d}}^{k-1}}{{\mathrm{d}}\,{\varOmega }^{k-1}}\left\{ ({\varOmega }-{\varOmega }_{1,2})^k\frac{[\chi ^{(3)}({\varOmega })]^k}{{\varOmega }-\omega '}\right\} \end{aligned}$$
(30)

For the first order poles (\(k=1\)) expression (30) gives the well-known result

$$\begin{aligned} {\mathrm{Res}}_{{\varOmega }={\varOmega }_{1,2}}\frac{\chi ^{(3)}({\varOmega })}{{\varOmega }-\omega '}=\lim _{{\varOmega }\rightarrow {\varOmega }_{1,2}}\left\{ ({\varOmega }-{\varOmega }_{1,2})\frac{\chi ^{(3)}({\varOmega })}{{\varOmega }-\omega '}\right\} . \end{aligned}$$
(31)

It is noteworthy here that the above residues have been calculated only for the powers of the susceptibility \(\omega ^m[\chi ^{(3)}(\omega )]^k\) i.e. with \(m=0\). This provides the fastest convergence, but one can naturally expand the calculation also for any converging value of m.

Fig. 2
figure 2

The real a and imaginary b parts of the degenerate third-order susceptibility, and the real c and imaginary d parts of the second power of the degenerate third-order susceptibility. The solid lines are exact values obtained from the classical anharmonic oscillator model from Eq. (4). The retrieved values are calculated by the conventional K–K relations (dashed lines) and by the generalized K–K relations (dots). In all cases the present analysis overperforms the conventional K–K analysis. The parameters used in the figures are: \(B=100\), \(\omega _0=3\times 10^{15}\) 1/s, and \({\varGamma }=1\times 10^{15}\) 1/s

Full size image

In Fig. 2a, b the solid lines represent the real and imaginary parts, respectively, of the degenerate third-order susceptibility, which are obtained from the classical anharmonic oscillator model given by Eq. (4). The results of the conventional K–K relations are shown by dashed lines and a qualitative fit to both spectra was observed. However, near the resonances the conventional K–K relations diverge substantially from the true spectra. This divergence was also observed to increase with a lower damping parameter \({\varGamma }\) value. A major improvement was observed as the effect of the residue terms was taken into account i.e. the poles that are located in the upper half plane and use the generalized K–K relations (dots) given by Eqs. (22) and (23) with \(m=1\) and \(k=1\). Here an analytic expression for the residue terms of the first-order poles was calculated using Eqs. (4) and (31) and the generalized K–K relations give an excellent fit throughout the spectra independent of the used damping parameter \({\varGamma }\) value. In Fig. 2c, d the real and imaginary parts are presented, respectively, of the second power of the degenerate susceptibility (solid lines). The calculations using the second power of the degenerate susceptibility require simultaneous knowledge of both real and imaginary part of the susceptibility. The conventional K–K relations (dashed lines) now gave a much better agreement for the second power although the fit was not perfect. This is encouraging since it allows us to utilize the second and higher powers of the meromorphic susceptibility in testing the validity of measured data by the conventional K–K relations. However, it is noteworthy that both real and imaginary part of the complex function is needed in the second power calculation that can be achieved, for example, by using the maximum entropy method (MEM) [3]. Secondly, typically only the modulus of the susceptibility can be measured in nonlinear optical spectroscopy. In such case the real and imaginary parts of the meromorphic susceptibility can also be obtained using the MEM. Thus, the validity of the real and imaginary parts obtained by MEM can be, to some extent, tested by exploiting the conventional K–K relations for the powers of the meromorphic complex susceptibility.

An additional benefit of the generalized K–K relations is the rapid convergence of these dispersion integrals that can be evaluated within a finite frequency range. The better agreement between the conventional and generalized K–K relations in Fig. 2c, d was caused by a faster convergence of the second power of the susceptibility (proportional to \(\omega ^{-16}\)) and the effect of the residue terms was reduced. However, after a lengthy algebra an analytic expression for the residue terms was obtained for the second power from the anharmonic oscillator model and an excellent fit was observed between the generalized K–K relations (dots) and the exact curves (solid lines). The effect of the second-order poles were taken into account by Eq. (30) with \(k=2\). In nonlinear terahertz spectroscopy [36] it is possible to measure both amplitude and phase of a wave that interacts nonlinearly with the medium. Hence, in the case of degenerate susceptibility in THz spectral range it is possible to test the validity of measured real and imaginary parts of the degenerate susceptibility by utilizing the dispersion theory of this article. Actually, the power of modified K–K relations was demonstrated in testing complex refractive index obtained from experimental transmission spectra of active pharmaceutical ingredients (API) in terahertz time-domain spectroscopy [37, 38].

Finally, K–K type dispersion relations given by Eqs. (22)–(29) are utilized in the derivation of sum rules for the degenerate third-order susceptibility in the static limit by setting \(\omega '=0\). The sum rules for degenerate nonlinear susceptibility are given in analogy with (16) and their physical interpretation is considered below. For even integers m with \(m\le 8k-2\) one obtains from Eqs. (22) and (23)

$$\begin{aligned}&\int _0^{\infty }\omega ^{m-1}{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega )]^k\}{\mathrm{d}}\omega \nonumber \\&\quad =\pi {{\mathrm{Re}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\{{\varOmega }^{m-1}[\chi ^{(3)}({\varOmega })]^k\}\right\} , \end{aligned}$$
(32)
$$\begin{aligned} {{\mathrm{Im}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\{{\varOmega }^{m-1}[\chi ^{(3)}({\varOmega })]^k\}\right\} =0, \end{aligned}$$
(33)

and for odd integers m with \(m\le 8k-1\) sum rules from Eqs. (27) and (28) are of the form

$$\begin{aligned}&\int _0^{\infty }\omega ^{m-1}{{\mathrm{Re}}}\{[\chi ^{(3)}(\omega )]^k\}{\mathrm{d}}\omega \nonumber \\&\quad =-\pi {{\mathrm{Im}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\{{\varOmega }^{m-1}[\chi ^{(3)}({\varOmega })]^k\}\right\} , \end{aligned}$$
(34)
$$\begin{aligned} {{\mathrm{Re}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\{{\varOmega }^{m-1}[\chi ^{(3)}({\varOmega })]^k\}\right\} =0. \end{aligned}$$
(35)

For the highest convergent value \(m=8k\) the sum rule given by Eq. (33) remains in the same form, but the imaginary part sum rule Eq. (32) is transformed according to the high frequency value of a susceptibility and from Eq. (29) one obtains

$$\begin{aligned}&\int _0^{\infty }\omega ^{m-1}{{\mathrm{Im}}}\{[\chi ^{(3)}(\omega )]^k\}{\mathrm{d}}\omega \nonumber \\&\quad =\pi {{\mathrm{Re}}}\left\{ \sum _{{{\mathrm{Im}}}{\varOmega }>0}^{\mathrm{{poles}}}{\mathrm{Res}}\{{\varOmega }^{m-1}[\chi ^{(3)}({\varOmega })]^k\}\right\} -\frac{\pi }{2}\psi ^k. \end{aligned}$$
(36)

These sum rules (32)–(36) yield useful information about degenerate susceptibility in the static limit. The meromorphicity of the degenerate susceptibility produces residue terms to sum rules, which are now presented for the first time for the moments and the powers of the degenerate third-order susceptibility \({\varOmega }^m[\chi ^{(3)}({\varOmega })]^k\).

Consider first sum the rule of Eq. (32). If \(m=0\) and \(k=1\), the left hand-side integral resembles the expression on the right hand-side in Eq. (16). However, while the imaginary part of the complex refractive index namely \(k(\omega )\) in Eq. (16) is related to a holomorphic (analytic) complex function, the function \({{\mathrm{Im}}}\chi ^{(3)}(\omega )\}\) is related to a meromorphic function. From Eqs. (4) and (5) one has \(\chi ^{(3)}(0)={{\mathrm{Re}}}\chi ^{(3)}(0)=B[D(0)]^{-4}\}\), i.e. a value one would get if \(\chi ^{(3)}(\omega )\}\) was a holomorphic function. However, this is not the case but yet at least under the assumption of an anharmonic oscillator model the residue terms can be calculated with the aid of Eqs. (4), (5) and (32). Actually the resonance frequencies and the line widths can be estimated from the measured spectra. Then it is possible to find the numerical values of the complex poles of the degenerate susceptibility. Hence, using measured data and theoretical model for the susceptibility, a sum rule can be calculated like the one given in Eq. (32). In the case of the anharmonic model the sum rule yields

$$\begin{aligned} {\mathrm{P}}\int _0^{\infty }\frac{{{\mathrm{Re}}}\chi ^{(3)}(\omega )}{\omega }{\mathrm{d}}\omega =\pi {{\mathrm{Re}}}\left[ \frac{B}{{\varOmega }_1-{\varOmega }_2}\left( \frac{D^3({\varOmega }_1)}{{\varOmega }_1}-\frac{D^3({\varOmega }_2)}{{\varOmega }_2}\right) \right] . \end{aligned}$$
(37)

This type of sum rule is practical because it has clear physical meaning namely using measured \({{\mathrm{Re}}}\chi ^{(3)}(\omega )\) it is possible to get information on the parameter B that, as already mentioned, is related to the electron number density of the nonlinear medium. Note that in nonlinear terahertz spectroscopy of the degenerate four wave mixing there will be a nonlinear contribution to the refractive index \(n(0)-1\), which could be estimated using a sum rule like the one given by Eq. (37). Sum rule (33) has no direct physical interpretation because no spectra is involved and same holds for sum rule of Eq. (35), whereas sum rule (34) involves spectral data and can be treated as a pair sum rule of (37) in case that \(m=0\) and \(k=1\).

The results of this article complement previous studies of the generalized K–K relations and sum rules in nonlinear optical- and THz-spectroscopy [39]. A possible new application of the dispersion theory of degenerate four wave mixing susceptibility is related to the estimation of the van der Waals interaction forces between nanoparticles under the influence of a strong electromagnetic radiation and utilizing the nonlinear phenomenon of four wave mixing. A key factor for the van der Waals interaction forces is the knowledge of the Hamaker’s constant [40], which can be obtained from the imaginary part of the dielectric function of a system. Usually the Hamaker’s constant is approximated using optical constants, which are obtained at different spectral ranges such as UV, VIS, IR, and micro waves. An important parameter used in the calculation of the Hamaker’s constant is the permittivity at a zero frequency \(\varepsilon (0)\) [40]. The Hamaker’s constant and the f-sum rule have also a close connection [41]. Here it is important to point out that the interaction forces experienced by the nanoparticles will be different in the case of four wave mixing than in the case of linear phenomena. A rigorous calculation of the Hamaker’s constant involves inversion of the imaginary part of the permittivity but using a purely imaginary frequency i.e. in these notations \(\omega '\) is replaced by \({\mathrm{i}}\omega '\) in the dispersion relations (21)–(29). A nonlinear contribution is observed to the total permittivity \(\varepsilon =\varepsilon _{\mathrm{{L}}}+\varepsilon _{\mathrm{{NL}}}\) when dealing with a nonlinear optical phenomenon of four wave mixing, where ”L” denotes the linear and ”NL” the nonlinear part, respectively. The Hamaker’s constant should be calculated using the imaginary part of the total permittivity, if the nonlinear light interaction is properly taken into account. However, this is a topic of another investigation and beyond the scope of this paper.

5 Conclusions

In this theoretical paper the moments and the powers of the third-order degenerate susceptibility were investigated that is a meromorphic function in the complex frequency plane. Generalized K–K type dispersion relations were presented for both even and odd moments and have given the rigorous upper limits of the converging moments. The theorem of residues was utilized to take into account the order of the poles. Moreover, static limit sum rules were also derived that contain residue terms caused by the meromorphic nature of the investigated susceptibility.

The developed K–K type dispersion relations and sum rules have remarkably faster convergence than the previously presented counterparts in the literature. Therefore, these dispersion relations and sum rules are assumed to be significant especially in testing of theoretical models describing the degenerate third-order susceptibility. Furthermore, it was observed that the conventional K–K relations are more accurate for the second and higher powers of the meromorphic susceptibilities due to the faster convergence of the integrals and the effect of the poles is reduced. Hence, one can use the conventional K–K relations to estimate the real and imaginary parts of the higher powers of the meromorphic susceptibility but the most reliable results are obtained from the generalized K–K relations. It is a straightforward procedure to extend the K–K type dispersion relations and sum rules, which are now presented for the third-order degenerate susceptibility, to arbitrary-order meromorphic nonlinear susceptibilities. Unfortunately, the application of the developed theory to measured spectrum remains cumbersome since at resonance points of the spectrum the knowledge of both the real and imaginary parts of the degenerate nonlinear susceptibility would be required.

It is believed that these results will find applications in nonlinear terahertz spectroscopy. In THz spectroscopy, linear or nonlinear, it is possible to get information on both the phase and the amplitude of the wave. Therefore, the generalized dispersion relations of this study can be used in confirming the validity of experimentally obtained real and imaginary parts of the meromorphic THz susceptibility.