Approximate Calculation of Functional Integrals Generated by Nonrelativistic and Relativistic Hamiltonians

Abstract

The discussion revolves around the most recent outcomes in the realm of approximating functional integrals through calculations. Review of works devoted to the application of functional integrals in quantum mechanics and quantum field theory, nuclear physics and in other areas is presented. Methods obtained by the authors for approximate calculation of functional integrals generated by nonrelativistic Hamiltonians are given. One of the methods is based on the expansion in eigenfunctions of the Hamiltonian. In an alternate approach, the functional integrals are tackled using the semiclassical approximation. Methods for approximate evaluation of functional integrals generated by relativistic Hamiltonians are presented. These are the methods using functional polynomial approximation (analogue of formulas of a given degree of accuracy) and methods based on the expansion in eigenfunctions of the Hamiltonian, generating a functional integral.

1. Introduction

Functional integrals are widely used in quantum mechanics and quantum field theory nuclear physics and in many other areas (see [1,2,3,4,5,6,7,8,9,10,11] and references therein). In particular, functional integrals are used in the description of quantum mechanical systems of particles with different Hamiltonians.

The wide application of functional integrals stimulates the development of methods for their approximate calculation. Significant progress has been made over the past half century in the field of approximate calculation of functional integrals.

The motivation for writing this article was the desire to review well-known existing methods and methods obtained by the authors in recent years for calculating functional integrals.

There are various methods for computing different types of functional integrals. Since the “Feynman measure” is not a measure in the mathematically strict sense, many different approaches to Feynman integrals have arisen, substantiating their construction and offering appropriate methods for their approximate calculation [12,13,14,15].

Among the methods used, we note the Monte Carlo method for the approximate calculation of functional integrals [16,17,18].

In the realm of employed techniques, the Monte Carlo method stands out as a means to approximately compute functional integrals (as indicated in references [16,17,18]. This approach hinges on expressing the target integral as the mathematical expectation of a random variable. Its arithmetic mean of independent realizations yields an estimated value for the integral in question.

There are methods for the approximate calculation of functional integrals, which are based on the discretization of space and time [19,20,21]. The theory on the space-time lattice leads to the need to calculate integrals of large (≥${10}^5$) multiplicity, which requires the use of powerful computers with large RAM and high speed.

Significant progress has been made in the field of methods for the approximate calculation of functional integrals, in which discretization is not required, and the argument of the integrand remains an element of the original space. Within the framework of the approach that does not require discretization, a method has been developed for the approximate calculation of functional integrals based on the use of approximate formulas that are exact on the class of functional polynomials of a given degree (see [11,16,22,23,24,25] and bibliographic lists to them). Such formulas are called formulas of a given degree of accuracy and are widely used to approximate the calculation of functional integrals.

For integrals arising from Hamiltonians, an additional computational approach can be introduced in conjunction with the existing methods. This approach relies on utilizing the Feynman–Kac formula [2], which provides an integral representation for the evolution operator’s kernel. It involves expressing this kernel through the eigenvalues and eigenvectors of the Hamiltonian.

Note that the formulas for a given degree of accuracy are especially effective for functional integrals over the space of functions given on short intervals. For integrals determined over large time intervals, it is necessary to take high precision formulas, which lead to high computational costs. The method based on the eigenfunction expansion of the Hamiltonian is effective for functional integrals defined over large time intervals. For small time intervals, it is necessary to use a large number of terms in the expansion, which leads to high computational costs.

To calculate functional integrals with respect to the conditional Wiener measure, the semiclassical approximation is also considered [9]. The semiclassical approximation uses the expansion of the action with respect to the classical trajectory [1], which can be interpreted as an expansion in powers of Planck’s constant.

Section 2 is devoted to the calculation of functional integrals generated by nonrelativistic Hamiltonians: the Hamiltonian of the harmonic oscillator, the Hamiltonian of the anharmonic oscillator, the Hamiltonian of a one-dimensional rectangular potential well and the Hamiltonian containing the Davidson potential. Functional integrals generated by the nonrelativistic Hamiltonians are integrals with respect to the conditional Wiener measure. For such integrals, the theory of formulas of a given degree accuracy is developed in detail and there are good review papers. Therefore, in Section 2, emphasis is put on the methods obtained by the authors for calculating integrals based on the expansion in terms of eigenfunctions of the Hamiltonian [26,27] and the method of semiclassical approximation of functional integrals [28].

In Section 3, the calculation of functional integrals generated by relativistic Hamiltonians is analyzed. This section presents the results obtained by the methods of the authors based on approximation by polynomials (analogous to formulas of a given degree of accuracy) [29] and the methods based on expansion in terms of eigenfunctions of the Hamiltonian [30,31]. For the nonrelativistic Hamiltonian, the equation for eigenvalues and eigenfunctions has the form of the stationary Schrödinger equation for which an approximate solution can be found. For the relativistic Hamiltonian, in the general case, it is not possible to obtain an equation for the eigenvalues and eigenfunctions of the Hamiltonian in the form of the stationary Schrödinger equation. Therefore, the problem of the stationary Schrödinger equation for eigenvalues and eigenfunctions is treated in the presence of a small parameter $1/\left(m{c}^2\right)$. We also consider the construction of an equation of the Schrödinger type for the eigenvalues and eigenfunctions in the case of pseudospin symmetry (the sum of the vector and scalar potentials is equal to zero).

Numerical results and comparison of approximate and exact values are presented in both Section 2 and Section 3.

All the mentioned methods presented in Section 2 and Section 3 were proposed by the authors of this review and their detailed description can be found in the original papers [26,27,28,29,30,31].

2. Calculation of Functional Integrals Generated by Nonrelativistic Hamiltonians

2.1. Method Based on Expansion in Eigenfunctions of the Hamiltonian

In this section, we consider the problem of approximate calculation of functional integrals with respect to the conditional Wiener measure of the form

$$\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d{\mu }_{x_s,x_t}\left(x\right)=$$

$$=\underset{{\displaystyle \underset{1\le j\le n}{max}}\Delta t_j\to 0}{lim}{\int }_R\underset{\cdots }{(n-1)}{\int }_{R}exp\left\{-\sum _{j=1}^n(t_j-t_{j-1})V\left(x_j\right)\right\}\times$$

$$\prod _{j=1}^n{K}_{t_j-t_{j-1}}^0(x_{j-1},x_j)d{x}_1\dots d{x}_{n-1},$$

(1)

if this limit exists for any partition of the interval $[s,t]$ by points $s=t_0<t_1<\cdots <t_n=t$.

Here, $x_j=x\left(t_j\right);$ $K_{t_j-t_{j-1}}^0(x_{j-1},x_j)=\frac{1}{\sqrt{2\pi (t_j-t_{j-1})}}exp\{-\frac{{(x_j-x_{j-1})}^2}{2(t_j-t_{j-1})}\}$ is the kernel of the operator $exp\left\{t{H}_0\right\}$, $H_0=\frac{1}{2}\frac{{\partial }^2}{\partial x^2}$. The functional integral in Equation (1) is said to be generated by the Hamiltonian $H=\frac{1}{2}\frac{{\partial }^2}{\partial x^2}-V$.

In Equation (1), the argument of the exponential is not the action functional, but only part of it corresponding to the potential energy. There is no kinetic term ${\left(\frac{dx}{dt}\right)}^2$ under the integral sign since it is included in the measure of integration. The measure of integration on cylindrical sets $\left\{x\right(\tau ):$ $x\left(s\right)=x_s$, $x\left(t\right)=x_t$, $x\left(t_j\right)\in I_j$, $j=1$, 2, …, $n-1\}$, $I_j$ are Borel sets from R, using functions $K_{t_j-t_{j-1}}^0(x_{j-1},x_j)$ is defined by the expression

$${\int }_{I_{n-1}}\underset{\cdots }{(n-1)}{\int }_{I_1}\prod _{j=1}^n{K}_{t_j-t_{j-1}}^0(x_{j-1},x_j)d{x}_1\dots d{x}_{n-1}.$$

That is, the functional integral defined by Equation (1) is an integral over the countably additive conditional Wiener measure in a mathematically rigorous sense [2].

In this paper, to calculate functional integrals in Equation (1), we consider a new method. This method is based on the use of the Feynman–Kac formula [2], which gives an integral representation for the kernel $K(x_s,x_t)$ of the evolution operator $exp\left\{tH\right\}=exp\left\{t(H_0-V)\right\},$ and express the kernel $K(x_s,x_t)$ with the use of eigenvalues and eigenvectors of the Hamiltonian H.

$K(x_s,x_t)$ satisfies the equality [2]

$$K(x_s,x_t)=\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d{\mu }_{x_s,x_t}\left(x\right).$$

(2)

The kernel $K(x_s,x_t)$ can also be represented as [32]

$$K(x_s,x_t)=\sum _{n=0}^{\infty }{\psi }_n\left(x_s\right){\psi }_n\left(x_t\right)exp\{-{\lambda }_n(t-s)\},$$

where $-{\lambda }_n$, ${\psi }_n\left(x\right)$ are the eigenvalues and eigenvectors of the operator $H=\frac{1}{2}\frac{{\partial }^2}{\partial x^2}-V$.

Thus, the following statement is true.

Statement 1.

The functional integral defined by Equation (1) satisfies the equality

$$\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d{\mu }_{x_s,x_t}\left(x\right)=\sum _{n=0}^{\infty }{\psi }_n\left(x_s\right){\psi }_n\left(x_t\right)exp\{-{\lambda }_n(t-s)\}.$$

To calculate $-{\lambda }_n$, ${\psi }_n\left(x\right)$ we consider functions on the interval $[-A\le x\le A],$ A is some positive number, and the operator $H=\frac{1}{2}\frac{{\partial }^2}{\partial x^2}-V$ is approximated by a finite difference operator with matrix $\overline{H}$ with dimension $(N-1)\times (N-1)$, resulting from the approximation of the second derivative at the node $x_j$ by the expression $h^{-2}(x_{j+1}-2x_j+x_{j-1})$

$$\overline{H}=\left(\begin{array}{ccccc}{b}_1& c_1& \cdots & 0& 0\\ a_1& b_2& \cdots & 0& 0\\ ⋮& ⋮& \ddots & ⋮& ⋮\\ 0& 0& \cdots & b_{N-2}& c_{N-2}\\ 0& 0& \cdots & a_{N-2}& b_{N-1}\end{array}\right)=$$

$$=\frac{1}{h^2}\left(\begin{array}{ccccc}-1-h^2{V}_1& \frac{1}{2}& 0& \cdots & 0\\ \frac{1}{2}& -1-h^2{V}_2& \frac{1}{2}& \cdots & 0\\ 0& \frac{1}{2}& -1-h^2{V}_3& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& 0& \cdots & -1-h^2{V}_{N-1}\end{array}\right),$$

(3)

where $V_j=V\left(x_j\right)=V(-A+jh),1\le j\le N-1,h=\frac{2A}{N}.$

To calculate the eigenvalues $-{\overline{\lambda }}_n$ of a three-diagonal matrix $\overline{H},$ one can use the Sturm sequence method [33]. Consider the following sequence of polynomials, known as the Sturm sequence

$$D_0\left(\lambda \right)=1,D_1\left(\lambda \right)=b_1-\lambda ,D_j\left(\lambda \right)=(b_j-\lambda )D_{j-1}\left(\lambda \right)-a_{j-1}{c}_{j-1}{D}_{j-2}\left(\lambda \right),$$

$2\le j\le N-1$.

For an arbitrary $\lambda$, we define the function $s\left(\lambda \right)$ as the number of coincidences of signs in successive members of the sequence

$$D_0\left(\lambda \right),D_1\left(\lambda \right),D_2\left(\lambda \right),\dots ,D_{N-1}\left(\lambda \right),$$

Moreover, if $D_j\left(\lambda \right)=0$, we take the sign of $D_{j-1}\left(\lambda \right)$ as the sign of this term. Then the value of the function $s\left(\lambda \right)$ is equal to the number of eigenvalues of the matrix $\overline{H}$ greater than or equal to $\lambda$.

To calculate the eigenvectors ${\overline{\psi }}_j$ of the three-diagonal matrix $\overline{H}$, one can use the backward iteration method [33]. Using the computed $-{\lambda }_n$, ${\psi }_n\left(x\right)$, $0\le n\le m$, and the approximate formula

$$\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d{\mu }_{x_s,x_t}\left(x\right)\approx \sum _{n=0}^m{\psi }_n\left(x_s\right){\psi }_n\left(x_t\right)exp\{-{\lambda }_n(t-s)\}$$

(4)

we obtain the numerical values of the functional integral in Equation (1).

2.2. Harmonic Oscillator

In this section, we consider the approximate calculation of the functional integral in Equation (1) in the case $V\left(x\right)=\frac{1}{2}(-1+x^2)$. This integral is generated by the Hamiltonian of the harmonic oscillator

$$H=\frac{1}{2}\frac{{\partial }^2}{\partial x^2}-\frac{1}{2}(-1+x^2).$$

For a given potential $V\left(x\right)$, using the Sturm sequence method, we obtain approximate values for eigenvalues $-{\overline{\lambda }}_0$, $-{\overline{\lambda }}_1$ of three diagonal matrix $\overline{H}$. For $A=4$, $N=20$, $-{\overline{\lambda }}_0=0.0022$, $-{\overline{\lambda }}_1=-0.9736$. For comparison, exact values are $-{\lambda }_0=0$, $-{\lambda }_1=-1$.

Using the backward iteration method, we find approximations of the eigenvectors ${\overline{\psi }}_0,{\overline{\psi }}_1$ of three diagonal matrix $\overline{H}$.

Figure 1, Figure 2, Figure 3 and Figure 4 for $A=4$, $N=20$ show the exact and approximate values for the functions ${\psi }_0\left(x\right)$ and ${\psi }_1\left(x\right)$ together with the corresponding absolute and relative errors of calculations, respectively. The exact values for ${\psi }_0\left(x\right),{\psi }_1\left(x\right)$ are obtained from the expressions

$${\psi }_0\left(x\right)=\sqrt[4]{\frac{1}{\pi }}exp\left\{-\frac{x^2}{2}\right\},{\psi }_1\left(x\right)=-\sqrt[4]{\frac{1}{\pi }}\sqrt{2}xexp\left\{-\frac{x^2}{2}\right\}.$$

Using the calculated eigenvalues $-{\lambda }_0$, $-{\lambda }_1$, eigenvectors ${\psi }_0$, ${\psi }_1$, and Equation (4), we obtain the numerical values of the functional integral.

In Figure 5, Figure 6, Figure 7 and Figure 8 below for $A=4$, $N=20$, $s=0$ and different values $x_s=x_t$, exact and approximate functional integral (1) values are given at $t=1$ and $t=5$, together with the corresponding absolute and relative errors of calculations, respectively. The exact values for the functional integral at $s=0$ are obtained from the expression [2,32]

$$K(x_s,x_t)=\frac{exp\{t/2\}}{\sqrt{2\pi sinh\left(t\right)}}exp\left\{-\frac{exp\left\{t\right\}{(exp\{-t\}{x}_t-x_s)}^2}{2sinh\left(t\right)}-x_t^2+x_s^2\right\}.$$

Figure 5 and Figure 7 show that the use of two terms in Equation (4) gives a good approximation of the functional integral.

2.3. Anharmonic Oscillator

In this section, we consider the approximate calculation of the functional integral in Equation (1) in the case $V\left(x\right)=\frac{1}{2}(x^2+x^4)$. This integral is generated by the Hamiltonian of the anharmonic oscillator

$$H=\frac{1}{2}\frac{{\partial }^2}{\partial x^2}-\frac{1}{2}(x^2+x^4).$$

Using the calculated eigenvalues $-{\lambda }_0$, $-{\lambda }_1$, eigenvectors ${\psi }_0$, ${\psi }_1$ and Equation (4), we obtain the numerical values of the functional integral. In Figure 9, Figure 10, Figure 11 and Figure 12 for $A=4$, $N=40$, $s=0$, $m=0$ and $m=1$, approximate values of the functional integral in Equation (1) are given for various values of $x_s=x_t$ at $t=1$ and $t=2$, together with the corresponding absolute and relative errors of calculations, respectively.

2.4. One-Dimensional Rectangular Potential Well

In this section, we consider the approximate calculation of the functional integral in Equation (1) in the case $V\left(x\right)=-{\pi }^2/\left(8a\right)$ at $-a<x<a$, $V\left(x\right)=+\infty$ at $\left|x\right|\ge a$. This integral is generated by the Hamiltonian of a one-dimensional rectangular well with absolutely rigid walls

$$H=\frac{1}{2}\frac{{\partial }^2}{\partial x^2}-V\left(x\right).$$

Using the calculated eigenvalues $-{\lambda }_0$, $-{\lambda }_1$, eigenvectors ${\psi }_0$, ${\psi }_1$, and Equation (4), we obtain the numerical values of the functional integral. In Figure 13, Figure 14, Figure 15 and Figure 16 for a $a=1$, $N=40$, $s=0$, $m=0$, and $m=1$, approximate values of the functional integral in Equation (1) are given for various values of $x_s=x_t$ at $t=1$ and $t=2$, together with the corresponding absolute and relative errors of calculations, respectively.

2.5. Davidson Potential

When describing particles in quantum mechanics, potentials are used that contain a variable in the denominator, for example, the Coulomb potential $V\left(x\right)=\frac{-c}{x}$, the shifted Coulomb potential $V\left(x\right)=-c_1/x+c_2$, the Kratzer potential $V\left(x\right)=c_1/x^2-c_2/x+c_3$, the Davidson potential $V\left(x\right)=c_1{x}^2+c_2/x^2$ and others [34]. The expression $c_2/x^2$ can be interpreted as a repulsive potential due to the centrifugal force. Therefore, this expression is usually called the centrifugal potential [35].

In this section, we consider the approximate calculation of the functional integral in Equation (1) in the case $V\left(x\right)=\frac{1}{2}(\frac{5}{16x^2}+x^2)$, where $V\left(x\right)=+\infty$ for $x\le 0$. That is, the Davidson potential is considered with the coefficients $c_1=1/2$, $c_2=5/32$.

Using the calculated eigenvalues $-{\lambda }_0$, $-{\lambda }_1$, eigenvectors ${\psi }_0$, ${\psi }_1$, and Equation (4), we obtain the numerical values of the functional integral.

Figure 17 and Figure 18 show the exact and approximate values of the integral $K(x_s,x_t)$ defined by Equation (2), $x_s=x_t=x$ at $s=0$, $N=40$, $m=0$ and $m=1$ for $t=1$, together with the corresponding absolute and relative errors of calculations, respectively.

Figure 19, Figure 20, Figure 21, Figure 22, Figure 23 and Figure 24 show the exact and approximate values of the integral $K(x_s,x_t)$ defined by Equation (2), $x_s=x_t=x$ at $s=0$, $m=1$, $N=10$, $N=40$ for $t=1$ and $t=2$, together with the corresponding absolute and relative errors of calculations, respectively.

Figure 17 and Figure 21 demonstrate that for various values of t, the approximate values are closer to the exact values for larger m. Figure 19 and Figure 23 show that for various values of t, the approximate values are closer to the exact values for larger N.

The exact values for the functional integral are obtained from the formula [36,37,38]

$$K(x_s,x_t)=\int exp\left\{-{\int }_s^t(c_2{x}^2\left(\tau \right)+\frac{c_3}{x^2\left(\tau \right)})d\tau \right\}d{\mu }_{x_s,x_t}\left(x\right)=$$

$$=\frac{\gamma \sqrt{x_s{x}_t}}{2sinh\left(\frac{\gamma (t-s)}{2}\right)}exp\left\{-\frac{\gamma }{4}(x_s^2+x_t^2)coth(\frac{\gamma (t-s)}{2})\right\}{I}_{\mu }\left(\frac{\gamma x_s{x}_t}{2sinh\frac{\gamma (t-s)}{2}}\right),$$

where $\gamma =\sqrt{8c_2},$ $\mu =\frac{1}{2}\sqrt{1+8c_3}$, $I_{\mu }$ are the modified Bessel functions of order $\mu$.

2.6. Semiclassical Approximation for the Calculation of Functional Integrals

This section deals with functional integrals

$$\int exp\left\{-\frac{1}{\hslash }{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d{\mu }_{x_s,x_t}\left(x\right),$$

(5)

which are defined by (1), where

$$K_{t_j-t_{j-1}}^0(x_{j-1},x_j)=\frac{1}{\sqrt{2\pi \hslash (t_j-t_{j-1})}}exp\left\{-\frac{{(x_j-x_{j-1})}^2}{2\hslash (t_j-t_{j-1})}\right\}$$

is the kernel of $exp\left\{\frac{t}{\hslash }{H}_0\right\}$, $H_0=\frac{{\hslash }^2}{2}\frac{{\partial }^2}{\partial x^2}$, ℏ is the parameter that takes positive real values.

Using the Feynman–Kac formula [2], the kernel $K_{t-s}(x_s,x_t)$ of the operator evolution $exp\left\{\frac{t}{\hslash }H\right\}$, $H=\frac{{\hslash }^2}{2}\frac{{\partial }^2}{\partial x^2}-V\left(x\right)$, is represented through the functional integral

$$K_{t-s}(x_s,x_t)=\int exp\left\{-\frac{1}{\hslash }{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d{\mu }_{x_s,x_t}\left(x\right).$$

In this section, we consider a semiclassical approximation for calculating the indicated functional integrals [9,39,40]. The semiclassical approximation uses the expansion of the action with respect to the classical trajectory [1], which can be interpreted as an expansion in powers of the Planck constant. In this section, we describe the method of semiclassical approximation of functional integrals and perform a numerical analysis of the accuracy of the semiclassical approximation of functional integrals. For numerical analysis, we compare the results obtained with the help of semiclassical approximation and with the help of the method of calculating functional integrals based on the expansion in terms of eigenfunctions of the Hamiltonian generating the functional integral. When approximating the functional integrals based on the expansion in terms of the eigenfunctions of the Hamiltonian, one can estimate the accuracy of the obtained approximate values. The semiclassical approximation can be used for a fairly wide range of problems. For example, the semiclassical approximation is used in the calculation of functional integrals used in field theory, in particular, in the theory of gravitation. However, to estimate the accuracy of the semiclassical approximation, a comparison with the results obtained by other methods is necessary.

2.6.1. Description of the Method

To calculate the functional integral in Equation (5), rewrite it as

$$K_{t-s}(x_s,x_t)=\int D\left[x\right]exp\left\{-\frac{1}{\hslash }{\int }_s^t\left[\frac{1}{2}{\left(\dot{x}\left(\tau \right)\right)}^2+V\left(x\left(\tau \right)\right)\right]d\tau \right\},$$

(6)

where

$$D\left[x\right]=\underset{n\to \infty }{lim}\prod _{i=1}^{n-1}\frac{d{x}_i}{\sqrt{2\pi \hslash (t_i-t_{i-1})}}\frac{1}{\sqrt{2\pi \hslash (t_n-t_{n-1})}}.$$

The quantity $\frac{1}{2}{\left(\dot{x}\left(\tau \right)\right)}^2+V\left(x\left(\tau \right)\right)$ can be considered as the Lagrangian of the system $L(\dot{x},x,\tau )$, and the value $S={\int }_s^{t}L(\dot{x},x,\tau )d\tau$ as the action. Using the principle of least action [1] it is possible to single out from all possible trajectories the classical trajectory $x_{\mathrm{cl}}$, for which the action S takes an extremal value. The classical trajectory is found as a solution of the Euler–Lagrange equations

$$\frac{d}{dt}\left[\frac{\partial L}{\partial \dot{x}}\right]-\frac{\partial L}{\partial x}=0.$$

Further, to calculate the integral, expansion of the action S with respect to the classical trajectory $x_{\mathrm{cl}}$ can be used

$$S\left[x\left(\tau \right)\right]\approx S\left[x_{\mathrm{cl}}\left(\tau \right)\right]+\frac{1}{2}{\delta }^{2}S\left[x_{\mathrm{cl}}\left(\tau \right)\right].$$

The second variation ${\delta }^{2}S\left[x_{\mathrm{cl}}\left(\tau \right)\right]$ can be written as

$${\delta }^{2}S\left[x_{\mathrm{cl}}\left(\tau \right)\right]={\int }_s^t\delta x\mathsf{\Lambda }\delta xd\tau ,$$

where $x=x_{\mathrm{cl}}+\delta x,$

$$\mathsf{\Lambda }=\frac{{\partial }^{2}L}{\partial x^2}{|}_{x=x_{\mathrm{cl}}}+\frac{{\partial }^{2}L}{\partial x\partial \dot{x}}{|}_{x=x_{\mathrm{cl}}}\frac{d}{dt}-\frac{d}{dt}\frac{{\partial }^{2}L}{\partial \dot{x}\partial x}{{|}_{x=x_{\mathrm{cl}}}-\frac{d}{dt}\frac{{\partial }^{2}L}{\partial {\dot{x}}^2}|}_{x=x_{\mathrm{cl}}}\frac{d}{dt}.$$

Thus, the integral in Equation (6) can be written in the form

$$K_{t-s}(x_s,x_t)=exp\left\{-\frac{1}{\hslash }S\left[x_{\mathrm{cl}}\left(\tau \right)\right]\right\}\int D\left[y\right]exp\left\{-\frac{1}{2\hslash }{\int }_s^{t}y\mathsf{\Lambda }yd\tau \right\},$$

(7)

where integration is performed over trajectories $y=\delta x$, satisfying the conditions $y\left(s\right)=0,$ $y\left(t\right)=0,$

$$D\left[y\right]=\underset{n\to \infty }{lim}\prod _{i=1}^{n-1}\frac{d{y}_i}{\sqrt{2\pi \hslash (t_i-t_{i-1})}}\frac{1}{\sqrt{2\pi \hslash (t_n-t_{n-1})}}.$$

To calculate the integral in Equation (7), we use the expansion

$$\sum _{j=1}^{\infty }{a}_j{u}_j,$$

where the functions $u_j$ are eigenfunctions of the operator $\mathsf{\Lambda }$ with eigenvalues ${\lambda }_j.$ Then, the integral in Equation (7) can be written as

$$K=\int D\left[y\right]exp\left\{-\frac{1}{2\hslash }{\int }_s^{t}y\mathsf{\Lambda }yd\tau \right\}=J\int D\left[a\right]exp\left\{-\frac{1}{2\hslash }\sum _{j=1}^{\infty }{\lambda }_j{a}_j^2\right\},$$

(8)

where J is the Jacobian of the transition from the variable y to the variable a,

$$D\left[a\right]=\underset{n\to \infty }{lim}\prod _{i=1}^{n-1}d{a}_i.$$

Since the Jacobian J is invariant under the operators $\mathsf{\Lambda }$ [1,40], then

$$K\prod _{j=1}^{\infty }{\lambda }_j^{\frac{1}{2}}=K_{\mathrm{free}}\prod _{j=1}^{\infty }{\lambda }_{\mathrm{free},j}^{\frac{1}{2}},$$

where

$$K_{\mathrm{free}}=J\int D\left[a\right]exp\left\{-\frac{1}{2\hslash }\sum _{j=1}^{\infty }{\lambda }_{\mathrm{free},j}{a}_j^2\right\},$$

and ${\lambda }_{\mathrm{free},j}$ are the eigenvalues of the operator $\mathsf{\Lambda }=-\frac{d^2}{d{t}^2}$.

Similar to Equation (8)

$$K_{\mathrm{free}}=\int D\left[y\right]exp\left\{-\frac{1}{2\hslash }{\int }_s^{t}y{\mathsf{\Lambda }}_{\mathrm{free}}yd\tau \right\}=$$

$$=\underset{n\to \infty }{lim}\int \cdots \int \prod _{i=1}^{n-1}\frac{d{y}_i}{\sqrt{2\pi \hslash (t_i-t_{i-1})}}\frac{1}{\sqrt{2\pi \hslash (t_n-t_{n-1})}}exp\left\{-\sum _{i=0}^{n-1}\frac{{(y_{i+1}-y_i)}^2}{2\hslash (t_{i+1}-t_i)}\right\}=$$

$$=\frac{1}{\sqrt{2\pi \hslash (t-s)}}.$$

Consequently,

$$K=K_{\mathrm{free}}\prod _{j=1}^{\infty }\frac{{\lambda }_{\mathrm{free},j}^{\frac{1}{2}}}{{\lambda }_j^{\frac{1}{2}}}=\frac{1}{\sqrt{2\pi \hslash (t-s)}}\prod _{j=1}^{\infty }\frac{{\lambda }_{\mathrm{free},j}^{\frac{1}{2}}}{{\lambda }_j^{\frac{1}{2}}}.$$

(9)

Thus, from Equations (8) and (9), it follows that the integral in Equation (7) is written in the form

$$K_{t-s}(x_s,x_t)=exp\left\{-\frac{1}{\hslash }S\left[x_{\mathrm{cl}}\left(\tau \right)\right]\right\}\frac{1}{\sqrt{2\pi \hslash (t-s)}}\prod _{j=1}^{\infty }\frac{{\lambda }_{\mathrm{free},j}^{\frac{1}{2}}}{{\lambda }_j^{\frac{1}{2}}}.$$

(10)

In our case $L(\dot{x},x,\tau )=\frac{1}{2}{\left(\dot{x}\left(\tau \right)\right)}^2+V\left(x\left(\tau \right)\right).$ Consequently, the trajectory $x_{\mathrm{cl}}$ satisfies the equation

$${\stackrel{..}{x}}_{\mathrm{cl}}\left(\tau \right)-V^{{}^{\prime }}\left(x_{\mathrm{cl}}\left(\tau \right)\right)=0.$$

(11)

Consequently,

$${\int }_s^t{\stackrel{.}{x}}_{\mathrm{cl}}^2\left(\tau \right)d\tau ={\stackrel{.}{x}}_{\mathrm{cl}}\left(t\right)x_{\mathrm{cl}}\left(t\right)-{\stackrel{.}{x}}_{\mathrm{cl}}\left(s\right)x_{\mathrm{cl}}\left(s\right)-{\int }_s^t{V}^{{}^{\prime }}\left(x_{\mathrm{cl}}\left(\tau \right)\right)x_{\mathrm{cl}}\left(\tau \right)d\tau ,$$

$$S\left[x_{\mathrm{cl}}\left(\tau \right)\right]={\int }_s^t\left[V\left(x_{\mathrm{cl}}\left(\tau \right)\right)-\frac{1}{2}{V}^{{}^{\prime }}\left(x_{\mathrm{cl}}\left(\tau \right)\right)x_{\mathrm{cl}}\left(\tau \right)\right]d\tau +$$

$$+\frac{1}{2}\left({\stackrel{.}{x}}_{\mathrm{cl}}\left(t\right)x_{\mathrm{cl}}\left(t\right)-{\stackrel{.}{x}}_{\mathrm{cl}}\left(s\right)x_{\mathrm{cl}}\left(s\right)\right).$$

(12)

Approximate values of the function $x_{\mathrm{cl}}\left(\tau \right),$ $s\le \tau \le t,$ are found by solving Equation (11) using the grid method for solving nonlinear boundary value problems [41].

The operator $\mathsf{\Lambda }$ in our case has the form

$$\mathsf{\Lambda }=V^{″}\left(x_{\mathrm{cl}}\left(\tau \right)\right)-\frac{d^2}{d{t}^2}.$$

This operator is approximated by a finite difference operator with a matrix $\overline{\mathsf{\Lambda }}$ of dimension $(N-1)\times (N-1),$ obtained as a result of approximating the second derivative at the node $t_j$ by the expression $\Delta t^{-2}(t_{j+1}-2t_j+t_{j-1}),$

$$\overline{\mathsf{\Lambda }}=\frac{1}{\Delta t^2}\left(\begin{array}{ccccc}2+\Delta t^2{V}_1& -1& 0& \cdots & 0\\ -1& 2+\Delta t^2{V}_2& -1& \cdots & 0\\ 0& -1& 2+\Delta t^2{V}_3& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& 0& \cdots & 2+\Delta t^2{V}_{N-1}\end{array}\right),$$

(13)

where $V_j=V^{″}\left(x_{\mathrm{cl}}(j\Delta t)\right),$ $1\le j\le N-1,$ $\Delta t=(t-s)/N$.

We replace the operator ${\mathsf{\Lambda }}_{\mathrm{free}}=-\frac{d^2}{d{t}^2}$ with a matrix of dimension $(N-1)\times (N-1)$

$${\overline{\mathsf{\Lambda }}}_{\mathrm{free}}=\frac{1}{\Delta t^2}\left(\begin{array}{ccccc}2& -1& 0& \cdots & 0\\ -1& 2& -1& \cdots & 0\\ 0& -1& 2& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& 0& \cdots & 2\end{array}\right).$$

Then, ${\prod }_{j=1}^{N-1}{\lambda }_j=det\overline{\mathsf{\Lambda }},$ ${\prod }_{j=1}^{N-1}{\lambda }_{\mathrm{free},j}=det{\overline{\mathsf{\Lambda }}}_{\mathrm{free}}.$

Thus, substituting approximate values for $S\left[x_{\mathrm{cl}}\left(\tau \right)\right]$ and ${\prod }_{j=1}^{N-1}{\lambda }_{\mathrm{free},j}^{\frac{1}{2}}/{\lambda }_j^{\frac{1}{2}}$ in Equation (10), we obtain a semiclassical approximation of the functional integral in Equation (5).

2.6.2. Numerical Results

In this section, we consider a numerical analysis of the accuracy of the semiclassical approximation of functional integrals using the example of an anharmonic oscillator, i.e., for the Hamiltonian $H=\frac{{\hslash }^2}{2}\frac{{\partial }^2}{\partial x^2}-V\left(x\right)$ and potential $V\left(x\right)=\frac{1}{2}{x}^2+\frac{1}{2}{x}^4.$ The functional integral has the form

$$K_{t-s}(x_s,x_t)=\int exp\left\{-\frac{1}{2\hslash }{\int }_s^t\left[x^2\left(\tau \right)+x^4\left(\tau \right)\right]d\tau \right\}d{\mu }_{x_s,x_t}\left(x\right).$$

(14)

In this case, $K_{t-s}(x_s,x_t)$ is calculated by Equation (10), where the trajectory $x_{\mathrm{cl}}$ satisfies the equation

$${\stackrel{..}{x}}_{\mathrm{cl}}\left(\tau \right)-x_{\mathrm{cl}}\left(\tau \right)-2x_{\mathrm{cl}}^3\left(\tau \right)=0.$$

From Equation (12), we obtain that

$$S\left[x_{\mathrm{cl}}\left(\tau \right)\right]=-\frac{1}{2}{\int }_s^t{x}_{\mathrm{cl}}^4\left(\tau \right)d\tau +\frac{1}{2}\left({\stackrel{.}{x}}_{\mathrm{cl}}\left(t\right)x_{\mathrm{cl}}\left(t\right)-{\stackrel{.}{x}}_{\mathrm{cl}}\left(s\right)x_{\mathrm{cl}}\left(s\right)\right).$$

Values $V_j$ in Equation (13) have the form $V_j=1+6x_{\mathrm{cl}}^2(j\Delta t),$ $1\le j\le N-1,$ $\Delta t=(t-s)/N$.

To estimate the error of the semiclassical approximation of the integral in Equation (5), we calculate this integral using an approximate method based on the expansion in terms of the eigenfunctions of the Hamiltonian. In this method, $K_{t-s}(x_s,x_t)$ is calculated by Equation (4), where ${\lambda }_j$ and ${\psi }_j$—eigenvalues and eigenvectors of the matrix in Equation (3), $V_j=\frac{1}{2}{(-A+jh)}^2+\frac{1}{2}{(-A+jh)}^4,$ $1\le j\le N-1.$

For comparison, Figure 25 shows the values of integral in Equation (14) at $s=0,$ $t=1,$ $\hslash =1,$ $x_s=0,$ $-3\le x_t\le 3,$ obtained using the semiclassical approximation and using the expansion in terms of the eigenfunctions of the Hamiltonian. Figure 26 shows the values of integral in Equation (14) at $s=0,$ $t=1,$ $\hslash =0.1,$ $x_s=0,$ $-1.5\le x_t\le 1.5,$ obtained using the semiclassical approximation and using the expansion in terms of eigenfunctions of the Hamiltonian.

Figure 25 and Figure 26 show that the semiclassical approximation approximates the functional integral well, and not only for small values $\hslash .$

Thus, in this section, a numerical analysis of the accuracy of the semiclassical approximation of functional integrals is carried out. For numerical analysis, we used a comparison of the results obtained using the semiclassical approximation and using the method of calculating functional integrals based on the expansion in terms of eigenfunctions of the Hamiltonian generating the functional integral.

3. Calculation of Functional Integrals Generated by Relativistic Hamiltonians

3.1. Method Based on Approximation by Polynomials

In this section, we consider matrix-valued functional integrals generated by the Hamiltonian $H=a\alpha \frac{\partial }{\partial x}+b\beta .$ Exact formulas for functional integrals of functional polynomials are proposed and methods based on approximation by polynomials are considered (an analogue of formulas with a given degree of accuracy). The results of a numerical experiment are also presented.

3.1.1. Analytical Results

The matrix-valued functional integral over the space of functions $x\left(\tau \right),\tau \in [0,t]$ is defined by the equality

$$\int F\left(x\right(·\left)\right)d\mu \left(x\right)=$$

(15)

$$=\underset{\underset{j}{max}\Delta t_j\to 0}{lim}\underset{{\mathbf{R}}^N}{\int }F\left(\sum _{j=1}^N{x}_j{\chi }_{]t_{j-1},t_j]}(·)\right)\prod _{j=N}^{1}S(\Delta t_j,\Delta x_j)d{x}_1\dots d{x}_N,$$

if this limit exists as $\underset{j}{max}\Delta t_j\to 0$ for any partition of the interval $[0,t]$ by points $0=t_0<t_1<\dots <t_N=t.$ Here, $\Delta t_j=t_j-t_{j-1},$ $\Delta x_j=x_j-x_{j-1},$ $x_j=x\left(t_j\right),j=1,\dots ,N,$ ${\chi }_{]t_{j-1},t_j]}(·)$ is the characteristic function of the interval $[t_{j-1},t_j].$ $S(t,x)$ is the fundamental solution of the equation [42,43]

$$\frac{\partial S(t,x)}{\partial t}=a\alpha \frac{\partial S(t,x)}{\partial x}+b\beta S(t,x),$$

(16)

where $a,b$ are real parameters, $\alpha ,\beta$ are anticommuting quantities (operators or matrices), i.e., $\alpha \beta +\beta \alpha =0.$ We also assume that ${\alpha }^2={\beta }^2=E$.

The following theorem is true for the integral (16) [29].

Theorem 1.

Let the functional integral be defined by Equation (15) and the transition function $S(t,x)$ be defined by Equation (16). Let also the functions $f_j\left(\tau \right),1\le j\le d$ be Riemann integrable on $[s,t].$ Then, the following equality holds:

$$\int {\int }_s^t{f}_1\left(\tau \right)dx\left(\tau \right)\cdots {\int }_s^t{f}_d\left(\tau \right)dx\left(\tau \right)d\mu \left(x\right)=$$

(17)

$$={(-a)}^d{\int }_s^t\left(d\right){\int }_s^t\prod _{k=1}^d{f}_k\left({\tau }_k\right)\prod _{k=0}^{d}exp\left\{({\tau }_{m_k}-{\tau }_{m_{k+1}}){(-1)}^{k}b\beta \right\}d{\tau }_1\dots d{\tau }_d{\alpha }^d,$$

where $(m_1,\dots ,m_d)$ is a permutation of $(1,\dots ,d)$ such that ${\tau }_{m_1}\ge {\tau }_{m_2}\ge \dots \ge {\tau }_{m_d},$ ${\tau }_{m_0}=t,$ ${\tau }_{m_{d+1}}=s.$

Let us give one more exact formula, which will be needed when comparing approximate values with exact ones [29].

Statement 2.

Let the functional integral be defined by Equation (15) and the transition function $S(t,x)$ be defined by Equation (16). Then, the following equality holds:

$$\int exp\left\{{\int }_s^{t}f\left(\tau \right)dx\left(\tau \right)\right\}d\mu \left(x\right)=\overleftarrow{T}exp\left\{{\int }_s^t\left(-af\left(\tau \right)\alpha +b\beta \right)d\tau \right\},$$

(18)

where the symbol $\overleftarrow{T}exp$ denotes the chronologically ordered exponential, which is defined by the equality

$$\overleftarrow{T}exp\left\{{\int }_s^t\left(-af\left(\tau \right)\alpha +b\beta \right)d\tau \right\}=\underset{\underset{j}{max}\Delta t_j\to 0}{lim}\prod _{j=N}^{1}exp\left\{\Delta t_j\left(-af\left({\tau }_{j-1}\right)\alpha +b\beta \right)\right\}.$$

Equations (17) and (18) are exact formulas for functional integrals. Equation (17) can be used for approximate calculation of the integrals

$$\int F\left({\int }_s^t{f}_1\left(\tau \right)dx\left(\tau \right),\dots ,{\int }_s^t{f}_d\left(\tau \right)dx\left(\tau \right)\right)d\mu \left(x\right).$$

Statement 3.

Let the functions $f_j\left(\tau \right),1\le j\le d$ be Riemann integrable on $[s,t].$ Let the function $F(u_1,\dots ,u_d)$ be expanded into a series ${\sum }_{p_1,\dots ,p_d=0}^{\infty }{F}_{p_1,\dots ,p_d}{u}_1^{p_1}\cdots u_d^{p_d}$ and the series ${\sum }_{p_1,\dots ,p_d=0}^{\infty }|F_{p_1,\dots ,p_d}\left|\right|u_1{|}^{p_1}\cdots {\left|u_d\right|}^{p_d}$ converges for $u_j=a{c}_j,$ $c_j>{\int }_s^t\left|f_j\left(\tau \right)\right|d\tau ,1\le j\le d$. Then,

$$\int F\left({\int }_s^t{f}_1\left(\tau \right)dx\left(\tau \right),\dots ,{\int }_s^t{f}_d\left(\tau \right)dx\left(\tau \right)\right)d\mu \left(x\right)=\sum _{p_1,\dots ,p_d=0}^{\infty }{F}_{p_1,\dots ,p_d}\times$$

(19)

$${(-a)}^P{\int }_s^t\left(P\right){\int }_s^t\prod _{k=1}^P{g}_k\left({\tau }_k\right)\prod _{k=0}^{P}exp\left\{({\tau }_{m_k}-{\tau }_{m_{k+1}}){(-1)}^{k}b\beta \right\}d{\tau }_1\dots d{\tau }_P{\alpha }^P$$

and the series converges in this equality. Here, $P=p_1+\dots +p_d,$ $g_i\left(\tau \right)=f_k\left(\tau \right),$ ${\sum }_{m=0}^{k-1}{p}_m\le i<{\sum }_{m=0}^k{p}_m,$ $p_0=1,$ $1\le k\le d,$ $(m_1,\dots ,m_P)$ is a permutation of$(1,\dots ,P)$ such that ${\tau }_{m_1}\ge {\tau }_{m_2}\ge \dots \ge {\tau }_{m_P},$ ${\tau }_{m_0}=t,$ ${\tau }_{m_{P+1}}=s.$ Proof follows from the expansion of the function $F(u_1,\dots ,u_d)$ and Equation (17). The following statement is also true [29].

Statement 4.

Let the functional integral be defined by Equation (15) and the transition function $S(t,x)$ be defined by Equation (16). Let $\alpha =\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right),$ $\beta =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right).$ Then, the elements of the integral $\int exp\left\{-{\int }_s^t{x}^2\left(\tau \right)d\tau \right\}d\mu \left(x\right)$ with indices $(1,2)$ and $(2,1)$ are equal to zero.

In the functional integral defined in Equation (15), the function $S(t,x)$ takes matrix values and the functional F takes real values. One can consider a functional integral for which both the function $S(t,x)$ and the functional F take matrix values. For example, the functional F has the form $exp\left\{{\int }_s^t\left[\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right]d\tau \right\}$, where $\beta =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right).$ In this case, the functional integral is defined in the following way:

$$\int exp\left\{{\int }_s^t\left[\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right]d\tau \right\}d\mu \left(x\right)=$$

$$=\underset{\underset{j}{max}\Delta t_j\to 0}{lim}{\int }_R\underset{\cdots }{\left(n\right)}{\int }_R\prod _{j=n}^1[exp\left\{(t_j-t_{j-1})\left[\beta \mathsf{\Phi }\left(x_j\right)+V\left(x_j\right)\right]\right\}\times$$

$$S(t_j-t_{j-1},x_j-x_{j-1})]d{x}_1\dots d{x}_n,$$

(20)

if this limit exists for any partition of the interval $[s,t]$ by points $s=t_0<t_1<\cdots <t_n=t$. Here, $x_j=x\left(t_j\right)$ and the transition function $S(t_j-t_{j-1},x_j-x_{j-1})$ satisfied Equation (16).

For these integrals, one can prove a statement similar to Statement 4.

Statement 5.

Let the functional integral be defined by (20) and the transition function $S(t,x)$ be defined by Equation (16), where $\alpha =\left(\begin{array}{cc}0& -1\\ 1& 0\end{array}\right),$ $\beta =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right).$ Let $\mathsf{\Phi }\left(x\right)=x^2/2$, $V\left(x\right)=-x^2/2$. Then, the elements of the integral

$$I=\int exp\left\{{\int }_s^t[\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)]d\tau \right\}d\mu \left(x\right)$$

with indices $(1,2)$ and $(2,1)$ are equal to zero.

Statements 4 and 5 will be needed to compare exact and approximate values.

3.1.2. Numerical Results

Before using Equation (19) for approximate calculations, consider the estimate for

$${(-a)}^P{\int }_s^t\left(P\right){\int }_s^t\prod _{k=1}^P{g}_k\left({\tau }_k\right)\prod _{k=0}^{P}exp\left\{\left({\tau }_{m_k}-{\tau }_{m_{k+1}}\right){(-1)}^{k}b\beta \right\}d{\tau }_1\dots d{\tau }_P.$$

Let the functions $|f_j\left(\tau \right)|\le C,1\le j\le d$ and $\beta =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right).$ Under the assumptions made, to evaluate the indicated expression, it suffices to estimate

$${\left|a\right|}^P{C}^P\left|{\int }_s^t\left(P\right){\int }_s^t\prod _{k=0}^{P}exp\left\{\left({\tau }_{m_k}-{\tau }_{m_{k+1}}\right){(-1)}^k\left|b\right|\right\}d{\tau }_1\dots d{\tau }_P\right|$$

and

$${\left|a\right|}^P{C}^P\left|{\int }_s^t\left(P\right){\int }_s^t\prod _{k=0}^{P}exp\left\{\left({\tau }_{m_k}-{\tau }_{m_{k+1}}\right){(-1)}^k(-|b\left|\right)\right\}d{\tau }_1\dots d{\tau }_P\right|.$$

For even P

$$\left|{\int }_s^t\left(P\right){\int }_s^t\prod _{k=0}^{P}exp\left\{\left({\tau }_{m_k}-{\tau }_{m_{k+1}}\right){(-1)}^k\left|b\right|\right\}d{\tau }_1\dots d{\tau }_P\right|=$$

$$=P!\left|{\int }_s^t{\int }_s^{{\tau }_P}\dots {\int }_s^{{\tau }_3}{\int }_s^{{\tau }_2}exp\left\{\left(t-2{\tau }_P+2{\tau }_{P-1}-\cdots -2{\tau }_2+2{\tau }_1-s\right)\left|b\right|\right\}d{\tau }_1\cdots d{\tau }_P\right|=$$

$$=\frac{P!}{2\left|b\right|}|{\int }_s^t{\int }_s^{{\tau }_P}\cdots {\int }_s^{{\tau }_4}{\int }_s^{{\tau }_3}[exp\left\{\right(t-2{\tau }_P+2{\tau }_{P-1}-\cdots +$$

$$+2{\tau }_3-s\left)\left|b\right|\right\}-exp\left\{\left(t-2{\tau }_P+2{\tau }_{P-1}-\cdots +-2{\tau }_2+s\right)\left|b\right|\right\}]d{\tau }_2\dots d{\tau }_P|\le$$

$$\frac{(t-s)P!}{2\left|b\right|}|{\int }_s^t{\int }_s^{{\tau }_P}\cdots {\int }_s^{{\tau }_5}{\int }_s^{{\tau }_4}exp\left\{\left(t-2{\tau }_P+2{\tau }_{P-1}-\cdots +2{\tau }_3-s\right)\left|b\right|\right\}d{\tau }_3\cdots d{\tau }_P|.$$

Continuing the estimates, we obtain as a result

$$\left|{\int }_s^t\left(P\right){\int }_s^t\prod _{k=0}^{P}exp\left\{\left({\tau }_{m_k}-{\tau }_{m_{k+1}}\right){(-1)}^k\left|b\right|\right\}d{\tau }_1\cdots d{\tau }_P\right|\le P!\frac{{(t-s)}^{P/2}}{{\left(2\right|b\left|\right)}^{P/2}}exp\left\{(t-s)\right|b\left|\right\}.$$

For odd P, a similar inequality is true

$$\left|{\int }_s^t\left(P\right){\int }_s^t\prod _{k=0}^{P}exp\{({\tau }_{m_k}-{\tau }_{m_{k+1}}){(-1)}^k\left|b\right|\}d{\tau }_1\dots d{\tau }_P\right|\le P!\frac{{(t-s)}^{\frac{P-1}{2}}}{{\left(2\right|b\left|\right)}^{\frac{P+1}{2}}}exp\left\{(t-s)\right|b\left|\right\}.$$

Estimates are obtained similarly: for even P

$$\left|{\int }_s^t\left(P\right){\int }_s^t\prod _{k=0}^{P}exp\left\{\left({\tau }_{m_k}-{\tau }_{m_{k+1}}\right){(-1)}^k(-|b\left|\right)\right\}d{\tau }_1\cdots d{\tau }_P\right|\le P!\frac{{(t-s)}^{\frac{P-2}{2}}}{{\left(2\right|b\left|\right)}^{\frac{P+2}{2}}}exp\left\{(t-s)\right|b\left|\right\},$$

for odd P

$$\left|{\int }_s^t\left(P\right){\int }_s^t\prod _{k=0}^{P}exp\left\{\left({\tau }_{m_k}-{\tau }_{m_{k+1}}\right){(-1)}^k(-|b\left|\right)\right\}d{\tau }_1\cdots d{\tau }_P\right|\le P!\frac{{(t-s)}^{\frac{P-1}{2}}}{{\left(2\right|b\left|\right)}^{\frac{P+1}{2}}}exp\left\{(t-s)\right|b\left|\right\}.$$

Since in the given estimates $\left|a\right|$ is in the numerator, and $\left|b\right|$ in the denominator, we find that Equation (19) is effective for small $\left|a\right|$ and large $\left|b\right|$.

As an example, consider the approximation of the integral

$$\int exp\left\{{\int }_s^t\lambda \left(\tau \right)dx\left(\tau \right)\right\}d\mu \left(x\right).$$

Using the expansion in Equation (19), we obtain that

$$\int exp\left\{{\int }_s^t\lambda \left(\tau \right)dx\left(\tau \right)\right\}d\mu \left(x\right)\approx \sum _{a_1=0}^n\frac{1}{a_1!}{(-a)}^{a_1}\times$$

$${\int }_s^t\left(a_1\right){\int }_s^t\prod _{k=1}^{a_1}\lambda \left({\tau }_k\right)\prod _{k=0}^{a_1}exp\left\{\left({\tau }_{m_k}-{\tau }_{m_{k+1}}\right){(-1)}^{k}b\beta \right\}d{\tau }_1\cdots d{\tau }_{a_1}{\alpha }^{a_1}=$$

$$=\sum _{a_1=0}^n{(-a)}^{a_1}{\int }_s^t\left(a_1\right){\int }_{{\tau }_3}^t{\int }_{{\tau }_2}^t\prod _{k=1}^{a_1}\lambda \left({\tau }_k\right)\times$$

$$exp\left\{\left(t-2{\tau }_1+2{\tau }_2-\cdots +{(-1)}^{a_1}2{\tau }_{a_1}-{(-1)}^{a_1}s\right)b\beta \right\}d{\tau }_1\cdots d{\tau }_{a_1}{\alpha }^{a_1}.$$

It follows from Equation (18) that for $\lambda \equiv 1$ the exact value of the integral is equal to

$$\int exp\left\{{\int }_s^t\lambda \left(\tau \right)dx\left(\tau \right)\right\}d\mu \left(x\right)=exp\left\{(t-s)(-a\alpha +b\beta )\right\}=$$

$$=cosh\left((t-s)\sqrt{a^2+b^2}\right)+\frac{-a\alpha +b\beta }{\sqrt{a^2+b^2}}sinh\left((t-s)\sqrt{a^2+b^2}\right).$$

For $t=2,$ $s=1,$ $\lambda \equiv 1$, $\alpha =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right),$ $\beta =\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right),$ $b=1$ and small $a=0.1$, the exact value of the integral is $\left(\begin{array}{cc}1.431& 1.177\\ 1.177& 1.667\end{array}\right).$

Approximate values for $n=0$ and $n=1$ are $\left(\begin{array}{cc}1.543& 1.175\\ 1.175& 1.543\end{array}\right),\left(\begin{array}{cc}1.426& 1.175\\ 1.175& 1.661\end{array}\right),$ respectively.

For $t=2,$ $s=1,$ $\lambda \equiv 1,$ $\alpha =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right),$ $\beta =\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right),$ $a=1$ and large $b=10$, the exact value of the integral is $\left(\begin{array}{cc}10414& 11507\\ 11507& 12715\end{array}\right).$

Approximate values for $n=0$ and $n=1$ are $\left(\begin{array}{cc}11014& 11014\\ 11014& 11014\end{array}\right),\left(\begin{array}{cc}9912& 11014\\ 11014& 12115\end{array}\right),$ respectively.

3.2. Functional Integrals Generated by the Matrix Equation

The method for calculating functional integrals based on the expansion in terms of eigenfunctions of the Hamiltonian generating the functional integral is applicable to the nonrelativistic Hamiltonian, since the equation for the eigenvalues and eigenfunctions of the nonrelativistic Hamiltonian has the form of the stationary Schrödinger equation, for which one can find an approximate solution. For the relativistic Hamiltonian, in the general case, it is not possible to obtain an equation for the eigenvalues and eigenfunctions of the Hamiltonian in the form of the stationary Schrödinger equation.

In this section, the derivation of the stationary Schrödinger equations for eigenvalues and eigenfunctions in the presence of a small parameter $1/\left(m{c}^2\right)$ will be considered. To obtain the stationary Schrödinger equation, an approximation is used, which is obtained by expanding into a series in powers of $1/\left(m{c}^2\right)$ up to terms of the smallest order.

3.2.1. Analytical Results

In this section, we consider matrix-valued functional integrals, which, following the works [42,43], are defined on the space of functions $x\left(\tau \right),$ $s\le \tau \le t,$ satisfying the condition $x\left(s\right)=0$ and the Lipschitz condition with order equal to one, that is for any $s\le r_1<r_2\le t,$ $|x\left(r_1\right)-x\left(r_2\right)|\le M|r_2-r_1.|$ The integrals are defined by the equality

$$\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)=\underset{\underset{j}{max}\Delta t_j\to 0}{lim}{\int }_R\underset{\cdots }{\left(n\right)}{\int }_{R}exp\left\{-\sum _{j=1}^n\left(t_j-t_{j-1}\right)V\left(x_j\right)\right\}\times$$

$$\prod _{j=n}^{1}S\left(t_j-t_{j-1},x_j-x_{j-1}\right)d{x}_1\cdots d{x}_n,$$

(21)

if this limit exists for any partition of the interval $[s,t]$ by points $s=t_0<t_1<\cdots <t_n=t.$

Here, $x_j=x\left(t_j\right);$ $S\left(t_j-t_{j-1},x_j-x_{j-1}\right)$ is the transition function, which is the fundamental solution of the equation [42,43]

$$\frac{\partial S(t,x)}{\partial t}=a\alpha \frac{\partial S(t,x)}{\partial x}+b\beta S(t,x),$$

(22)

where a, b are real parameters, $\alpha$, $\beta$ are anticommuting matrices, i.e., $\alpha \beta +\beta \alpha =0$.

This integral is related to the kernel of the evolution operator [42,43], namely

$$\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)={\int }_{-\infty }^{+\infty }K(0,x_t)d{x}_t,$$

where $K(x_s,x_t)$ is the kernel of the operator $exp\left\{tH\right\},$ $H=a\alpha \frac{\partial }{\partial x}+b\beta -V\left(x\right).$

Our goal is to use the expansion of the function $K(x_s,x_t)$ in terms of eigenfunctions of the Hamiltonian H to calculate the functional integral.

For the relativistic Hamiltonian $H,$ the parameters a and b have the form $a=-c$, $b=-m{c}^2$, where m is the particle mass and c is the speed of light [44]. Therefore, we consider the calculation of eigenvalues and eigenvectors of the operator H for large values $\left|b\right|.$ We will also assume that

$$\alpha =\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right),\beta =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right).$$

Under the above assumptions, the equations for the eigenvalues E and eigenfunctions ${\left(\phi \left(x\right),\psi \left(x\right)\right)}^T$ (T is the transposition sign) of the operator H have the form

$$a{\partial }_x\psi \left(x\right)+b\phi \left(x\right)-V\left(x\right)\phi \left(x\right)=E\phi \left(x\right),$$

$$a{\partial }_x\phi \left(x\right)-b\psi \left(x\right)-V\left(x\right)\psi \left(x\right)=E\psi \left(x\right).$$

Denote $E^{{}^{\prime }}=E+b$. Then,

$$\begin{array}{c}{\displaystyle \left(E^{{}^{\prime }}-2b+V\left(x\right)\right)\phi \left(x\right)=a{\partial }_x\psi \left(x\right),}\\ {\displaystyle \left(E^{{}^{\prime }}+V\left(x\right)\right)\psi \left(x\right)=a{\partial }_x\phi \left(x\right).}\end{array}$$

(23)

From these equalities, we obtain an equation for $\psi \left(x\right)$

$$\left(E^{{}^{\prime }}+V\left(x\right)\right)\psi \left(x\right)=a{\partial }_x{\left(E^{{}^{\prime }}-2b+V\left(x\right)\right)}^{-1}a{\partial }_x\psi \left(x\right).$$

(24)

For large $\left|b\right|$, the approximate equality is true [44]

$${\left(E^{{}^{\prime }}-2b+V\left(x\right)\right)}^{-1}\approx {(-2b)}^{-1}\left(1-\frac{E^{{}^{\prime }}+V\left(x\right)}{-2b}\right).$$

With the help of this equality, Equation (24) is transformed into

$$\left(E^{{}^{\prime }}+V\left(x\right)\right)\psi \left(x\right)=\frac{a^2}{-2b}\left(\left[1-\frac{E^{{}^{\prime }}+V\left(x\right)}{-2b}\right]{\partial }_x^2+\frac{V^{{}^{\prime }}\left(x\right)}{2b}{\partial }_x\right)\psi \left(x\right).$$

Multiplying both sides of this equality by ${\left[1-\frac{E^{{}^{\prime }}+V\left(x\right)}{-2b}\right]}^{-1}\approx 1-\left(E^{{}^{\prime }}+V\left(x\right)\right)/\left(2b\right),$ we obtain

$$\left[\frac{a^2}{-2b}{\partial }_x^2-\left(1+\frac{E^{{}^{\prime }}+V\left(x\right)}{-2b}\right)\left(E^{{}^{\prime }}+V\left(x\right)\right)+\frac{a^2}{-2b}\left(1+\frac{E^{{}^{\prime }}+V\left(x\right)}{-2b}\right)\frac{V^{{}^{\prime }}\left(x\right)}{2b}{\partial }_x\right]\psi \left(x\right)\approx$$

$$\approx \left[\frac{a^2}{-2b}{\partial }_x^2-\left(E^{{}^{\prime }}+V\left(x\right)\right)-\frac{{\left(E^{{}^{\prime }}+V\left(x\right)\right)}^2}{-2b}+\frac{a^2}{-2b}\frac{V^{{}^{\prime }}\left(x\right)}{2b}{\partial }_x\right]\psi \left(x\right)\approx 0.$$

(25)

To find the eigenvalues $E^{{}^{\prime }}$ and eigenfunctions $\psi \left(x\right)$, consider the unperturbed operator

$$\frac{a^2}{-2b}{\partial }_x^2-V\left(x\right)$$

and assume that the eigenvalues $E_{0n}^{{}^{\prime }}$ and eigenfunctions ${\psi }_{0n}\left(x\right)$ are known for it, that is, the exact solutions of the equation

$$\left(\frac{a^2}{-2b}{\partial }_x^2-V\left(x\right)\right){\psi }_{0n}\left(x\right)=E_{0n}^{{}^{\prime }}{\psi }_{0n}\left(x\right)$$

are known.

The eigenvalues $E_{0n}^{{}^{\prime }}$ will be calculated approximately using the Sturm sequence method [33]. The eigenfunctions ${\psi }_{0n}\left(x\right)$ will be calculated approximately by the backward iteration method [33]. Then we seek the eigenvalues $E_n^{{}^{\prime }}$ and eigenfunctions ${\psi }_n\left(x\right)$ of the perturbed operator with $ϵ=1/(-2b)$ in the form [45]

$$E_n^{{}^{\prime }}=E_{0n}^{{}^{\prime }}+ϵE_{1n}^{{}^{\prime }}+{ϵ}^2{E}_{2n}^{{}^{\prime }}+\cdots$$

$${\psi }_n\left(x\right)=\sum _{j=0}^{\infty }{c}_j{\psi }_{0j}\left(x\right),c_j=c_{0j}+ϵc_{1j}+{ϵ}^2{c}_{2j}+\cdots$$

(26)

For an eigenvector with number n, $c_{0j}=1$ for $j=n$, $c_{0j}=0$ for $j\ne n$.

To find the first order corrections to the eigenvalues and eigenvectors, we substitute Equation (26) into (25). As a result, we obtain

$$\sum _{j=0}\left(E_{0j}^{{}^{\prime }}-E_{0n}^{{}^{\prime }}-ϵE_{1n}^{{}^{\prime }}\right)\left(c_{0j}+ϵc_{1j}\right){\psi }_{0j}\left(x\right)-$$

$$-ϵ{\left(E_{0n}^{{}^{\prime }}+V\left(x\right)\right)}^2\sum _{j=0}{c}_{0j}{\psi }_{0j}\left(x\right)-ϵ\frac{a^2}{-2b}{V}^{{}^{\prime }}\left(x\right){\partial }_x\sum _{j=0}{c}_{0j}{\psi }_{0j}\left(x\right)=0.$$

Taking into account that

$$\sum _{j=0}\left(E_{0j}^{{}^{\prime }}-E_{0n}^{{}^{\prime }}\right)c_{0j}{\psi }_{0j}\left(x\right)=0$$

the last equality will be written in the form

$$\sum _{j=0}\left(E_{0j}^{{}^{\prime }}-E_{0n}^{{}^{\prime }}\right)c_{1j}{\psi }_{0j}\left(x\right)-E_{1n}^{{}^{\prime }}{\psi }_{0n}\left(x\right)-$$

$$-{\left(E_{0n}^{{}^{\prime }}+V\left(x\right)\right)}^2{\psi }_{0n}\left(x\right)-\frac{a^2}{-2b}{V}^{{}^{\prime }}\left(x\right){\partial }_x{\psi }_{0n}\left(x\right)=0.$$

(27)

We introduce the notation

$$V_{jn}=-{\int }_{-\infty }^{+\infty }{\psi }_{0j}\left(x\right){\left(E_{0n}^{{}^{\prime }}+V\left(x\right)\right)}^2{\psi }_{0n}\left(x\right)dx-\frac{a^2}{-2b}{\int }_{-\infty }^{+\infty }{\psi }_{0j}\left(x\right)V^{{}^{\prime }}\left(x\right){\partial }_x{\psi }_{0n}\left(x\right)dx.$$

The scalar product of vectors ${\psi }_{ij}\left(x\right),{\psi }_{kl}\left(x\right)$ is defined through the integral

$${\int }_{-\infty }^{+\infty }{\psi }_{ij}\left(x\right){\psi }_{kl}\left(x\right)dx.$$

Multiplying both sides of Equation (27) scalarly by the vector ${\psi }_{0n}\left(x\right)$, we obtain the expression for $E_{1n}^{{}^{\prime }}$

$$E_{1n}^{{}^{\prime }}=V_{nn}.$$

Multiplying both sides of Equation (27) scalarly by the vector ${\psi }_{0j}\left(x\right)$, we obtain the expression for $c_{1j}$, $j\ne n$

$$c_{1j}=\frac{V_{jn}}{E_{0n}^{{}^{\prime }}-E_{0j}^{{}^{\prime }}}.$$

The coefficient $c_{1n}$ must be chosen so that the function ${\psi }_{0n}\left(x\right)+ϵ{\sum }_{j=0}{c}_{1j}{\psi }_{0j}\left(x\right)$ is normalized up to and including terms of the first order. To this end it is necessary to put $c_{1n}=0$ [45].

Thus, with the first-order correction, the eigenvector is given by

$${\psi }_{1n}\left(x\right)={\psi }_{0n}\left(x\right)+ϵ\sum _{j\ne n}\frac{V_{jn}}{E_{0n}^{{}^{\prime }}-E_{0j}^{{}^{\prime }}}{\psi }_{0j}\left(x\right).$$

(28)

It follows from Equation (23) that the functions ${\phi }_{ij}\left(x\right)$, $i=0$, 1, $j=0$, 1, 2, … are expressed in terms of ${\psi }_{ij}\left(x\right)$ by the formula

$${\phi }_{ij}\left(x\right)=\frac{a}{-2b}\left[1-\frac{E_{ij}^{{}^{\prime }}+V\left(x\right)}{-2b}\right]{\partial }_x{\psi }_{ij}\left(x\right).$$

Thus, the functions ${({\phi }_{ij}\left(x\right),{\psi }_{ij}\left(x\right))}^T$ are eigenfunctions of the operator H with eigenvalues $E_{ij}^{{}^{\prime }}-b,$ where $i=0$, 1 are the orders of smallness, $j=0$, 1, 2, … are the numbers of eigenvalues.

For a fixed i, the vectors $({\phi }_{ij}\left(x\right)$, ${\psi }_{ij}{\left(x\right))}^T$ are generally not orthogonal for different values of j. For their orthogonalization, the Gram-Schmidt process can be used. In this case, orthogonal vectors are constructed according to the rule

$${\eta }_{i0}={({\phi }_{i0}\left(x\right),{\psi }_{i0}\left(x\right))}^T,{\eta }_{ij}={({\phi }_{ij}\left(x\right),{\psi }_{ij}\left(x\right))}^T-\sum _{k=0}^{j-1}\frac{\langle {({\phi }_{ij}\left(x\right),{\psi }_{ij}\left(x\right))}^T,{\eta }_{ik}\rangle }{\langle {\eta }_{ik},{\eta }_{ik}\rangle }{\eta }_{ik},$$

(29)

where $\langle {({\eta }_1\left(x\right),{\eta }_2\left(x\right))}^T,{({\xi }_1\left(x\right),{\xi }_2\left(x\right))}^T\rangle ={\int }_{-\infty }^{+\infty }{\eta }_1\left(x\right){\xi }_1\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\eta }_2\left(x\right){\xi }_2\left(x\right)dx.$

To expand the kernel of the operator $exp\left\{tH\right\}$ with respect to the vector ${\eta }_{ij}$, one can use the formula

$$K(x_s,x_t)=\sum _{j=0}{\langle {\eta }_{ij},{\eta }_{ij}\rangle }^{-1}exp\left\{tH\right\}{\eta }_{ij}\left(x_t\right){\eta }_{ij}^T\left(x_s\right).$$

Since, according to Equation (29), ${\eta }_{ij}$ can be written as ${\eta }_{ij}={\sum }_{k=0}^j{c}_k{({\phi }_{ik}\left(x\right),{\psi }_{ik}\left(x\right))}^T$, where $c_k$ are some constants, $exp\left\{tH\right\}{\eta }_{ij}\left(x_t\right)$ can be written as

$$exp\left\{tH\right\}{\eta }_{ij}\left(x_t\right)=\sum _{k=0}^j{c}_{k}exp\left\{t(E_{ik}^{{}^{\prime }}-b)\right\}{\left({\phi }_{ik}\left(x\right),{\psi }_{ik}\left(x\right)\right)}^T.$$

Thus, the following statement is true.

Statement 6.

The functional integral defined in Equation (21) satisfies the approximate formula

$$\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)\approx$$

(30)

$${\int }_{-\infty }^{+\infty }\sum _{j=0}{\langle {\eta }_{ij},{\eta }_{ij}\rangle }^{-1}\sum _{k=0}^j{c}_{k}exp\left\{t(E_{ik}^{{}^{\prime }}-b)\right\}{\left({\phi }_{ik}\left(x_t\right),{\psi }_{ik}\left(x_t\right)\right)}^T{\eta }_{ij}^T\left(0\right)d{x}_t,$$

where $E_{ij}^{{}^{\prime }}-b$ are the approximate eigenvalues of the Hamiltonian H, ${({\phi }_{ij}\left(x\right),{\psi }_{ij}\left(x\right))}^T$ are non-normalized approximate eigenfunctions of the Hamiltonian H, ${\eta }_{ij}$ are normalized approximate eigenfunctions of the Hamiltonian H, $i=0$, 1 are the orders of smallness, $j=0$, 1, 2, … are the numbers of eigenvalues, $c_k$ are the expansion coefficients of the normalized functions by non-normalized functions.

3.2.2. Numerical Results

As an example, consider the calculation of the functional integral in Equation (21) in the case when the transition function S is determined by Equation (22), $\alpha =\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right),$ $\beta =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right),$ $V\left(x\right)=x^2/8$.

For the zero order of smallness $(i=0)$ and zero eigenvalue $(j=0)$, Equation (30) has the form

$$\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)\approx {\int }_{-\infty }^{+\infty }{〈({\phi }_{00},{\psi }_{00}),({\phi }_{00},{\psi }_{00})〉}^{-1}\times$$

$$exp\left\{t(E_{00}^{{}^{\prime }}-b)\right\}{\left({\phi }_{00}\left(x_t\right),{\psi }_{00}\left(x_t\right)\right)}^T\left({\phi }_{00}\left(0\right),{\psi }_{00}\left(0\right)\right)d{x}_t=$$

$$={〈\left({\phi }_{00},{\psi }_{00}\right),\left({\phi }_{00},{\psi }_{00}\right)〉}^{-1}exp\left\{t(E_{00}^{{}^{\prime }}-b)\right\}\times$$

$$\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{00}\left(x_t\right)d{x}_t{\phi }_{00}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{00}\left(x_t\right)d{x}_t{\psi }_{00}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{00}\left(x_t\right)d{x}_t{\phi }_{00}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{00}\left(x_t\right)d{x}_t{\psi }_{00}\left(0\right)\end{array}\right).$$

For $s=0$, $t=2$, $a=-4$, $b=-16$, $m=1$, $A=5$ ($[-A,A]$ is the interval on which approximate eigenfunctions are considered), $N=100$ (N is the number of intervals into which the interval $[-A,A]$ is divided for approximate calculation of eigenvalues and eigenfunctions) $E_{00}^{{}^{\prime }}\approx -0.2499$, ${\int }_{-\infty }^{+\infty }{\psi }_{00}\left(x_t\right)d{x}_t\approx 2.2366$, ${\int }_{-\infty }^{+\infty }{\phi }_{00}{x}_t)d{x}_t\approx 0.00048$, ${\psi }_{00}\left(0\right)\approx 0.6318$, ${\phi }_{00}\left(0\right)\approx 0.002$, ${〈({\phi }_{00},{\psi }_{00}),({\phi }_{00},{\psi }_{00})〉}^{-1}\approx 0.9923$.

Thus, for the zero order of smallness $(i=0)$ and zero eigenvalue $(j=0)$

$$exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)\approx \left(\begin{array}{cc}{10}^{-6}& 2·{10}^{-4}\\ 2·{10}^{-3}& 0.8507\end{array}\right).$$

For the zero order of smallness $(i=0)$ and the first eigenvalue $(j=1)$, the term with $j=1$ in Equation (30) has the form

$${〈({\phi }_{01},{\psi }_{01}),({\phi }_{01},{\psi }_{01})〉}^{-1}exp\left\{t(E_{01}^{{}^{\prime }}-b)\right\}\times$$

$$\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{01}\left(x_t\right)d{x}_t{\phi }_{01}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{01}\left(x_t\right)d{x}_t{\psi }_{01}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{01}\left(x_t\right)d{x}_t{\phi }_{01}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{01}\left(x_t\right)d{x}_t{\psi }_{01}\left(0\right)\end{array}\right).$$

For $s=0$, $t=2$, $a=-4$, $b=-16$, $m=1$, $A=5$, $N=100$ $E_{01}^{{}^{\prime }}\approx -0.7498$, ${\int }_{-\infty }^{+\infty }{\psi }_{01}\left(x_t\right)d{x}_t\approx {10}^{-6}$, ${\int }_{-\infty }^{+\infty }{\phi }_{01}\left(x_t\right)d{x}_t\approx 0.01747$, ${\psi }_{01}\left(0\right)\approx 8·{10}^{-7}$, ${\phi }_{01}\left(0\right)\approx 0.08$, ${〈({\phi }_{01},{\psi }_{01}),({\phi }_{01},{\psi }_{01})〉}^{-1}\approx 0.9804$.

Thus, for the zero order of smallness $(i=0)$, the term with $j=1$ in Equation (30) for $exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)$ gives the contribution

$$\left(\begin{array}{cc}2·{10}^{-4}& 2·{10}^{-9}\\ 2·{10}^{-8}& 2·{10}^{-13}\end{array}\right).$$

For the zero order of smallness $(i=0)$ and the second eigenvalue $(j=1)$ at $s=0$, $t=2$, $a=-4$, $b=-16$, $m=1$, $A=5$, $N=100$, ${\eta }_{02}=\left({\eta }_{02,1},{\eta }_{02,2}\right)\approx {\left({\phi }_{02}\left(x\right),{\psi }_{02}\left(x\right)\right)}^T+0.0109{\left({\phi }_{00}\left(x\right),{\psi }_{00}\left(x\right)\right)}^T$, ${\langle {\eta }_{02},{\eta }_{02}\rangle }^{-1}\approx 0.9626$. The term with $j=2$ in Equation (30) has the form

$${\langle {\eta }_{02},{\eta }_{02}\rangle }^{-1}exp\left\{t(E_{02}^{{}^{\prime }}-b)\right\}\times \left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{02}\left(x_t\right)d{x}_t{\eta }_{02,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{02}\left(x_t\right)d{x}_t{\eta }_{02,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{02}\left(x_t\right)d{x}_t{\eta }_{02,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{02}\left(x_t\right)d{x}_t{\eta }_{02,2}\left(0\right)\end{array}\right)+$$

$$+0.0109{\langle {\eta }_{02},{\eta }_{02}\rangle }^{-1}exp\left\{t\left(E_{00}^{{}^{\prime }}-b\right)\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{00}\left(x_t\right)d{x}_t{\eta }_{02,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{00}\left(x_t\right)d{x}_t{\eta }_{02,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{00}\left(x_t\right)d{x}_t{\eta }_{02,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{00}\left(x_t\right)d{x}_t{\eta }_{02,2}\left(0\right)\end{array}\right).$$

$E_{02}^{{}^{\prime }}\approx -1.2506$, ${\int }_{-\infty }^{+\infty }{\psi }_{02}\left(x_t\right)d{x}_t\approx 1.5408$, ${\int }_{-\infty }^{+\infty }{\phi }_{02}\left(x_t\right)d{x}_t\approx 9·{10}^{-4}$, ${\psi }_{02}\left(0\right)\approx -0.4484$, ${\phi }_{02}\left(0\right)\approx 7·{10}^{-3}$, ${\eta }_{02,1}\left(0\right)\approx 7·{10}^{-3}$, ${\eta }_{02,2}\left(0\right)\approx -0.4421$, ${〈({\phi }_{02},{\psi }_{02}),({\phi }_{02},{\psi }_{02})〉}^{-1}$ $\approx 0.9626$.

Thus, for the zero order of smallness $(i=0)$, the term with $j=2$ in Equation (30) for $exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)$ gives the contribution

$$\left(\begin{array}{cc}5·{10}^{-7}& -3·{10}^{-5}\\ 9·{10}^{-4}& -0.0596\end{array}\right).$$

For the zero order of smallness $(i=0)$ and for the same values of the parameters, taking into account the contribution of three terms, we obtain

$$exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)\approx \left(\begin{array}{cc}2·{10}^{-4}& 2·{10}^{-4}\\ 2·{10}^{-3}& 0.7911\end{array}\right).$$

(31)

For the first order of smallness $(i=1)$ and zero eigenvalue $(j=0)$ Equation (30) has the form

$$\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)\approx {〈({\phi }_{10},{\psi }_{10}),({\phi }_{10},{\psi }_{10})〉}^{-1}\times$$

$$exp\left\{t\left(E_{10}^{{}^{\prime }}-b\right)\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{10}\left(x_t\right)d{x}_t{\phi }_{10}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{10}\left(x_t\right)d{x}_t{\psi }_{10}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{10}\left(x_t\right)d{x}_t{\phi }_{10}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{10}\left(x_t\right)d{x}_t{\psi }_{10}\left(0\right)\end{array}\right).$$

It follows from Equation (28) that for $s=0$, $t=2$, $a=-4$, $b=-16$, $m=1$, $A=5$, $N=100$

$${\psi }_{10}\left(x\right)\approx {\psi }_{00}\left(x\right)+{10}^{-4}{\psi }_{01}\left(x\right)+0.0014{\psi }_{02}\left(x\right),$$

$${\phi }_{10}\left(x\right)\approx {\phi }_{00}\left(x\right)+{10}^{-4}{\phi }_{01}\left(x\right)+0.0014{\phi }_{02}\left(x\right).$$

With the specified parameters

$$E_{10}^{{}^{\prime }}\approx -0.2494,{\int }_{-\infty }^{+\infty }{\psi }_{10}\left(x_t\right)d{x}_t\approx 2.2388,{\int }_{-\infty }^{+\infty }{\phi }_{10}\left(x_t\right)d{x}_t\approx 0.00048,$$

$${\psi }_{10}\left(0\right)\approx 0.6312,{\phi }_{10}\left(0\right)\approx 0.002,{〈({\phi }_{10},{\psi }_{10}),({\phi }_{10},{\psi }_{10})〉}^{-1}\approx 0.9923.$$

Thus, for the first order of smallness $(i=1)$ and zero eigenvalue $(j=0)$

$$exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)\approx \left(\begin{array}{cc}6·{10}^{-7}& 2·{10}^{-4}\\ 2·{10}^{-3}& 0.8516\end{array}\right).$$

For the first order of smallness $(i=1)$ and the first eigenvalue $(j=1)$ for the specified parameters, from Equation (28) expressing ${\psi }_{11}\left(x\right)$ in terms of ${\psi }_{0j}\left(x\right)$, $j=0$, 1, 2, we obtain

$${\psi }_{11}\left(x\right)\approx {\psi }_{01}\left(x\right)+3·{10}^{-4}{\psi }_{00}\left(x\right)-7·{10}^{-6}{\psi }_{02}\left(x\right),$$

$${\phi }_{11}\left(x\right)\approx {\phi }_{01}\left(x\right)+3·{10}^{-4}{\phi }_{00}\left(x\right)-7·{10}^{-6}{\phi }_{02}\left(x\right).$$

After orthogonalization, we obtain

$${\eta }_{11}=\left({\eta }_{11,1},{\eta }_{11,2}\right)\approx {\left({\phi }_{11}\left(x\right),{\psi }_{11}\left(x\right)\right)}^T-4·{10}^{-4}{\left({\phi }_{10}\left(x\right),{\psi }_{10}\left(x\right)\right)}^T,$$

${\langle {\eta }_{11},{\eta }_{11}\rangle }^{-1}\approx 0.9804$.

The term with $j=1$ in Equation (30) has the form

$${〈{\eta }_{11},{\eta }_{11}〉}^{-1}exp\left\{t\left(E_{11}^{{}^{\prime }}-b\right)\right\}\times \left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{11}\left(x_t\right)d{x}_t{\eta }_{11,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{11}\left(x_t\right)d{x}_t{\eta }_{11,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{11}\left(x_t\right)d{x}_t{\eta }_{11,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{11}\left(x_t\right)d{x}_t{\eta }_{11,2}\left(0\right)\end{array}\right)-$$

$$4·{10}^{-4}{\langle {\eta }_{11},{\eta }_{11}\rangle }^{-1}exp\left\{t\left(E_{10}^{{}^{\prime }}-b\right)\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{10}\left(x_t\right)d{x}_t{\eta }_{11,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{10}\left(x_t\right)d{x}_t{\eta }_{11,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{10}\left(x_t\right)d{x}_t{\eta }_{11,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{10}\left(x_t\right)d{x}_t{\eta }_{11,2}\left(0\right)\end{array}\right).$$

With the specified parameters

$$E_{11}^{{}^{\prime }}\approx -0.7552,{\int }_{-\infty }^{+\infty }{\psi }_{11}\left(x_t\right)d{x}_t\approx 6·{10}^{-4},{\int }_{-\infty }^{+\infty }{\phi }_{11}\left(x_t\right)d{x}_t\approx 0.01747,$$

$${\psi }_{11}\left(0\right)\approx 2·{10}^{-4},{\phi }_{11}\left(0\right)\approx 8·{10}^{-2},{\eta }_{11,1}\left(0\right)\approx 8·{10}^{-2},{\eta }_{11,2}\left(0\right)\approx -5·{10}^{-5}.$$

Thus, for the first order of smallness $(i=1)$, the term with $j=1$ in Equation (30) for $exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)$ gives the contribution

$$\left(\begin{array}{cc}3·{10}^{-4}& -2·{10}^{-7}\\ -3·{10}^{-5}& -3·{10}^{-8}\end{array}\right).$$

For the first order of smallness $(i=1)$ and the second eigenvalue (28) expressing ${\psi }_{12}\left(x\right)$ in terms of ${\psi }_{0j}\left(x\right)$, $j=0$, 1, 2, we obtain

$${\psi }_{12}\left(x\right)\approx {\psi }_{02}\left(x\right)-0.0069{\psi }_{00}\left(x\right)+0.0179{\psi }_{01}\left(x\right),$$

$${\phi }_{12}\left(x\right)\approx {\phi }_{02}\left(x\right)-0.0069{\phi }_{00}\left(x\right)+0.0179{\phi }_{01}\left(x\right).$$

After orthogonalization, we obtain

$${\eta }_{12}=({\eta }_{12,1},{\eta }_{12,2})\approx {({\phi }_{12}\left(x\right),{\psi }_{12}\left(x\right))}^T+0.016{({\phi }_{10}\left(x\right),{\psi }_{10}\left(x\right))}^T-0.018{({\phi }_{11}\left(x\right),{\psi }_{11}\left(x\right))}^T,$$

 ${\langle {\eta }_{12},{\eta }_{12}\rangle }^{-1}\approx 0.9626.$

The term with $j=2$ in Equation (30) has the form

$${\langle {\eta }_{12},{\eta }_{12}\rangle }^{-1}exp\left\{t\left(E_{12}^{{}^{\prime }}-b\right)\right\}\times \left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{12}\left(x_t\right)d{x}_t{\eta }_{12,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{12}\left(x_t\right)d{x}_t{\eta }_{12,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{12}\left(x_t\right)d{x}_t{\eta }_{12,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{12}\left(x_t\right)d{x}_t{\eta }_{12,2}\left(0\right)\end{array}\right)+$$

$$+0.016{\langle {\eta }_{12},{\eta }_{12}\rangle }^{-1}exp\left\{t\left(E_{10}^{{}^{\prime }}-b\right)\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{10}\left(x_t\right)d{x}_t{\eta }_{12,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{10}\left(x_t\right)d{x}_t{\eta }_{12,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{10}\left(x_t\right)d{x}_t{\eta }_{12,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{10}\left(x_t\right)d{x}_t{\eta }_{12,2}\left(0\right)\end{array}\right)-$$

$$-0.018{\langle {\eta }_{12},{\eta }_{12}\rangle }^{-1}exp\left\{t\left(E_{11}^{{}^{\prime }}-b\right)\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_{11}\left(x_t\right)d{x}_t{\eta }_{12,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_{11}\left(x_t\right)d{x}_t{\eta }_{12,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_{11}\left(x_t\right)d{x}_t{\eta }_{12,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_{11}\left(x_t\right)d{x}_t{\eta }_{12,2}\left(0\right)\end{array}\right).$$

With the specified parameters

$$E_{12}^{{}^{\prime }}\approx -1.2677,{\int }_{-\infty }^{+\infty }{\psi }_{12}\left(x_t\right)d{x}_t\approx 1.5254,{\int }_{-\infty }^{+\infty }{\phi }_{12}\left(x_t\right)d{x}_t\approx {10}^{-3},$$

$${\psi }_{12}\left(0\right)\approx -0.4528,{\phi }_{12}\left(0\right)\approx 8·{10}^{-3},{\eta }_{12,1}\left(0\right)\approx 7·{10}^{-3},{\eta }_{12,2}\left(0\right)\approx -0.4427.$$

Thus, for the first order of smallness $(i=1)$, the term with $j=2$ in Equation (30) for $exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)$ gives the contribution

$$\left(\begin{array}{cc}-2·{10}^{-7}& -2·{10}^{-6}\\ {10}^{-3}& -0.0607\end{array}\right).$$

For the first order of smallness $(i=1)$ and for the same values of the parameters, taking into account the contribution of three terms, we obtain

$$exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)\approx \left(\begin{array}{cc}2·{10}^{-4}& 2·{10}^{-4}\\ 3·{10}^{-3}& 0.7909\end{array}\right).$$

(32)

Thus, Equations (31) and (32) give approximate values of the functional integral for zero and first orders of smallness, respectively. These formulas take into account the first three eigenfunctions of the Hamiltonian H.

Statement 4 implies that for the elements of the matrix

$$exp\left\{tb\right\}\int exp\left\{-{\int }_s^{t}V\left(x\left(\tau \right)\right)d\tau \right\}d\mu \left(x\right)$$

with indices $(1,2)$ and $(2,1)$ for $V\left(x\right)=x^2/8$ the exact values are known and they are equal to zero. For comparison, Figure 27 and Figure 28 show the approximate values of the elements of the matrix-valued integral with indices $(1,2)$ and $(2,1)$ for various values N.

Figure 27 and Figure 28 show that the approximate values of the elements of the matrix-valued integral tend to the exact value as the number of intervals N increases into which the interval $[-A,A]$ is divided for the approximate calculation of eigenvalues and eigenfunctions.

3.3. Functional Integrals Generated by the Dirac Equation with Pseudospin Symmetry

As noted earlier, for the relativistic Hamiltonian, in the general case, it is not possible to obtain an equation for the eigenvalues and eigenfunctions of the Hamiltonian in the form of the stationary Schrödinger equation. In the previous section, we considered the derivation of the stationary Schrödinger equation for eigenvalues and eigenfunctions in the presence of a small parameter $1/\left(m{c}^2\right)$. To obtain the stationary Schrödinger equation an approximation was used, which is obtained by expanding in a series in powers of $1/\left(m{c}^2\right)$ up to terms of the smallest order.

In this section, we consider the case of pseudospin symmetry [46] (the sum of the vector and scalar potentials is equal to zero). In this case, an equation of the Schrödinger type is constructed for the eigenvalues and eigenfunctions. Then, a method is constructed for calculating functional integrals based on the expansion in terms of eigenfunctions of the relativistic Hamiltonian generating the functional integral.

3.3.1. Analytical Results

In this section, matrix-valued functional integrals are considered, which are defined on the space of functions $x\left(\tau \right)$, $s\le \tau \le t$, satisfying the condition $x\left(s\right)=0$. The integrals are defined by the equality

$$\int exp\left\{{\int }_s^t\left(\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right)d\tau \right\}d\mu \left(x\right)=$$

$$=\underset{\underset{j}{max}\Delta t_j\to 0}{lim}{\int }_R\underset{\cdots }{\left(n\right)}{\int }_R\prod _{j=n}^{1}exp\left\{\left(t_j-t_{j-1}\right)\left(\beta \mathsf{\Phi }\left(x_j\right)+V\left(x_j\right)\right)\right\}\times$$

$$\prod _{j=n}^{1}S\left(t_j-t_{j-1},x_j-x_{j-1}\right)d{x}_1\cdots d{x}_n,$$

(33)

if this limit exists for any partition of the interval $[s,t]$ by points $s=t_0<t_1<\cdots <t_n=t$.

Here, $x_j=x\left(t_j\right)$; $S\left(t_j-t_{j-1},x_j-x_{j-1}\right)$ is the transition function, which is the fundamental solution of the equation [42,43]

$$\frac{\partial S(t,x)}{\partial t}=a\alpha \frac{\partial S(t,x)}{\partial x}+b\beta S(t,x),$$

(34)

where a, b are real parameters and $\alpha$, $\beta$ are anticommuting matrices, i.e., $\alpha \beta +\beta \alpha =0$.

This integral can be written as [42,43]

$$\int exp\left\{{\int }_s^t\left(\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right)d\tau \right\}d\mu \left(x\right)={\int }_{-\infty }^{+\infty }K(0,x_t)d{x}_t,$$

(35)

where $K(x_s,x_t)$ is the kernel of the evolution operator $exp\left\{tH\right\}$ with the relativistic Hamiltonian $H=a\alpha \partial /\partial x+b\beta +\beta \mathsf{\Phi }\left(x\right)+V\left(x\right)$ corresponding to the Dirac equation.

The functional $\mathsf{\Phi }\left(x\right)$ is called the scalar potential, the functional $V\left(x\right)$ is called the vector potential. It follows from Equation (35) that to calculate the functional integral, one can use the expansion of the function $K(x_s,x_t)$ in terms of the eigenfunctions of the Hamiltonian H that generates the functional integral.

To satisfy the condition of pseudospin symmetry and anticommutability of matrices $\alpha$, $\beta$, we assume that

$$\alpha =\left(\begin{array}{cc}0& -1\\ 1& 0\end{array}\right),\beta =\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right).$$

Under the above assumptions, the equations for the eigenvalues E and eigenfunctions ${\left(\phi \left(x\right),\psi \left(x\right)\right)}^T$ (T is the transposition sign) of the operator H have the form

$$\begin{array}{c}-a{\partial }_x\psi \left(x\right)+b\phi \left(x\right)+\mathsf{\Phi }\left(x\right)\phi \left(x\right)+V\left(x\right)\phi \left(x\right)=E\phi \left(x\right),\\ a{\partial }_x\phi \left(x\right)-b\psi \left(x\right)-\mathsf{\Phi }\left(x\right)\psi \left(x\right)+V\left(x\right)\psi \left(x\right)=E\psi \left(x\right).\end{array}$$

(36)

When constructing an equation of the Schrödinger type for eigenvalues and eigenfunctions, we consider the case of pseudospin symmetry, i.e., $\mathsf{\Phi }\left(x\right)+V\left(x\right)=0$. Then, Equations (36) take the form

$$a{\partial }_x\psi \left(x\right)-\left(b-E\right)\phi \left(x\right)=0,$$

$$a{\partial }_x\phi \left(x\right)+\left(-b+2V\left(x\right)-E\right)\psi \left(x\right)=0.$$

From these equalities, we obtain the equation for $\psi \left(x\right)$

$$\left[{\partial }_x^2+\frac{(b-E)2V\left(x\right)}{a^2}-\frac{b^2-E^2}{a^2}\right]\psi \left(x\right)=0.$$

(37)

The function $\phi \left(x\right)$ is expressed in terms of the function $\psi \left(x\right)$ by the equality

$$\phi \left(x\right)=\frac{a}{b-E}{\partial }_x\psi \left(x\right).$$

(38)

To find $E_j$ and ${\psi }_j\left(x\right)$, satisfying Equation (37), we use the Sturm sequence method [33] to find an approximate eigenvalue ${\lambda }_j$ of the operator ${\partial }_x^2+\frac{(b-E)2V\left(x\right)}{a^2}$. Calculate $E_j^{{}^{\prime }}=-\sqrt{b^2-{\lambda }_j{a}^2}$. We select $E_j$ so that $E_j\approx E_j^{{}^{\prime }}$ with a given accuracy. These $E_j$ are chosen as approximate eigenvalues. Then, using the reverse iteration method [33], we find the eigenvectors ${\psi }_j\left(x\right)$.

Using Equation (38), we find functions ${\left({\phi }_j\left(x\right),{\psi }_j\left(x\right)\right)}^T$ that are eigenfunctions of the operator H with eigenvalues $E_j$.

The vectors ${\left({\phi }_j\left(x\right),{\psi }_j\left(x\right)\right)}^T$ are generally not orthogonal for different values j. To orthogonalize them, one can use the Gram–Schmidt process. In this case, orthogonal vectors are constructed according to the rule

$${\eta }_0={\left({\phi }_0\left(x\right),{\psi }_0\left(x\right)\right)}^T,{\eta }_j={\left({\phi }_j\left(x\right),{\psi }_j\left(x\right)\right)}^T-\sum _{k=0}^{j-1}\frac{〈{\left({\phi }_j\left(x\right),{\psi }_j\left(x\right)\right)}^T,{\eta }_k〉}{\langle {\eta }_k,{\eta }_k\rangle }{\eta }_k,$$

(39)

where $〈{\left({\eta }_1\left(x\right),{\eta }_2\left(x\right)\right)}^T,{\left({\xi }_1\left(x\right),{\xi }_2\left(x\right)\right)}^T〉={\int }_{-\infty }^{+\infty }{\eta }_1\left(x\right){\xi }_1\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\eta }_2\left(x\right){\xi }_2\left(x\right)dx$.

To expand the kernel of the operator $exp\left\{tH\right\}$ in vectors ${\eta }_j$, one can use the formula

$$K(x_s,x_t)=\sum _{j=0}{〈{\eta }_j,{\eta }_j〉}^{-1}\left[exp\left\{tH\right\}{\eta }_j\right]\left(x_t\right){\eta }_j^T\left(x_s\right).$$

From this relation and Equation (35) we obtain the expression for the functional integral

$$\begin{array}{c}{\displaystyle \int exp\left\{{\int }_s^t\left(\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right)d\tau \right\}d\mu \left(x\right)=}\\ {\displaystyle =\sum _{j=0}{\int }_{-\infty }^{+\infty }{〈{\eta }_j,{\eta }_j〉}^{-1}\left[exp\left\{tH\right\}{\eta }_j\right]\left(x_t\right){\eta }_j^T\left(0\right)d{x}_t.}\end{array}$$

(40)

Since ${\eta }_j$ are expressed in terms of the eigenvectors of the operator H, $exp\left\{tH\right\}{\eta }_j$ can be written in terms of $exp\left\{t{E}_k\right\}$ and ${\left({\phi }_k\left(x\right),{\psi }_k\left(x\right)\right)}^T$, $0\le k\le j$.

Thus, the following statement is true.

Statement 7.

The functional integral defined by Equations (33) and (34) satisfies Equation (40), where ${\left({\phi }_j\left(x\right),{\psi }_j\left(x\right)\right)}^T$ are unnormalized eigenfunctions of the Hamiltonian H, ${\eta }_j$ are normalized eigenfunctions of the Hamiltonian H.

3.3.2. Numerical Results

Let us consider the application of the proposed method to the calculation of the functional integral in Equation (35) where $\mathsf{\Phi }\left(x\right)=x^2/2$, $V\left(x\right)=-x^2/2$.

If we consider only the term with $j=0$ in Equation (40), we obtain

$$\int exp\left\{{\int }_s^t\left[\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right]d\tau \right\}d\mu \left(x\right)\approx$$

$$\approx {\int }_{-\infty }^{+\infty }{〈{\eta }_0,{\eta }_0〉}^{-1}\left[exp\left\{tH\right\}{\eta }_0\right]\left(x_t\right){\eta }_0^T\left(0\right)d{x}_t,$$

 ${\eta }_0={\left({\phi }_0\left(x\right),{\psi }_0\left(x\right)\right)}^T$. Therefore, this approximate equality can be rewritten as

$$\int exp\left\{{\int }_s^t\left[\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right]d\tau \right\}d\mu \left(x\right)\approx {\left({\int }_{-\infty }^{+\infty }{\phi }_0^2\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\psi }_0^2\left(x\right)dx\right)}^{-1}\times$$

$$exp\left\{t{E}_0\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_0\left(x_t\right)d{x}_t{\phi }_0\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_0\left(x_t\right)d{x}_t{\psi }_0\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_0\left(x_t\right)d{x}_t{\phi }_0\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_0\left(x_t\right)d{x}_t{\psi }_0\left(0\right)\end{array}\right).$$

For $s=0$, $t=1$, $a=1$, $b=1$, $A=5$ ($[-A,A]$ is the interval on which approximate eigenfunctions are considered), $N=90$ (N is the number of intervals into which the interval $[-A,A]$ is divided for approximate calculation of eigenvalues and eigenfunctions)

$$E_0\approx -1.617,{\int }_{-\infty }^{+\infty }{\psi }_0\left(x_t\right)d{x}_t\approx 1.669,{\int }_{-\infty }^{+\infty }{\phi }_0\left(x_t\right)d{x}_t\approx 0,$$

$${\psi }_0\left(0\right)\approx 0.8478,{\phi }_0\left(0\right)\approx -0.0292,{\int }_{-\infty }^{+\infty }{\phi }_0^2\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\psi }_0^2\left(x\right)dx\approx 1.1181.$$

Thus,

$$\int exp\left\{{\int }_s^t\left[\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right]d\tau \right\}d\mu \left(x\right)\approx \left(\begin{array}{cc}0& 0\\ -0.0087& 0.2512\end{array}\right).$$

Let us now calculate the contribution from the term with $j=1$ in Equation (40). In this case, the vector ${({\phi }_1\left(x\right),{\psi }_1\left(x\right))}^T$ is orthogonal to the vector ${({\phi }_0\left(x\right),{\psi }_0\left(x\right))}^T.$ Consequently, ${\eta }_1={({\phi }_1\left(x\right),{\psi }_1\left(x\right))}^T.$ Therefore, the term with $j=1$ in Formula (41) can be written as

$${\left[{\int }_{-\infty }^{+\infty }{\phi }_1^2\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\psi }_1^2\left(x\right)dx\right]}^{-1}\times$$

$$exp\left\{t{E}_1\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_1\left(x_t\right)d{x}_t{\phi }_1\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_1\left(x_t\right)d{x}_t{\psi }_1\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_1\left(x_t\right)d{x}_t{\phi }_1\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_1\left(x_t\right)d{x}_t{\psi }_1\left(0\right)\end{array}\right).$$

For $s=0$, $t=1$, $a=1$, $b=1$, $A=5$, and $N=90$

$$E_1\approx -2.581,{\int }_{-\infty }^{+\infty }{\psi }_1\left(x_t\right)d{x}_t\approx {10}^{-5},{\int }_{-\infty }^{+\infty }{\phi }_1\left(x_t\right)d{x}_t\approx 0,$$

$${\psi }_1\left(0\right)\approx 4·{10}^{-6},{\phi }_1\left(0\right)\approx -0.4753,{\int }_{-\infty }^{+\infty }{\phi }_1^2\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\psi }_1^2\left(x\right)dx\approx 1.5745.$$

Thus, the term with $j=1$ contributes

$$\left(\begin{array}{cc}0& 0\\ -2·{10}^{-7}& 2·{10}^{-12}\end{array}\right).$$

Let us now calculate the contribution from the term with $j=2$ in Equation (40). The vector ${\left({\phi }_2\left(x\right),{\psi }_2\left(x\right)\right)}^T$ is orthogonal to the vector ${\left({\phi }_1\left(x\right),{\psi }_1\left(x\right)\right)}^T$, but not orthogonal to the vector ${\left({\phi }_0\left(x\right),{\psi }_0\left(x\right)\right)}^T$. Therefore, based on Equation (39) ${\eta }_2=\left({\eta }_{2,1},{\eta }_{2,2}\right)\approx {\left({\phi }_2\left(x\right),{\psi }_2\left(x\right)\right)}^T+0.081{\left({\phi }_0\left(x\right),{\psi }_0\left(x\right)\right)}^T$. Consequently, the term with $j=2$ in Equation (40) can be written as

$${\left[{\int }_{-\infty }^{+\infty }{\eta }_{2,1}^2\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\eta }_{2,2}^2\left(x\right)dx\right]}^{-1}\times$$

$$exp\left\{t{E}_2\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_2\left(x_t\right)d{x}_t{\eta }_{2,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_2\left(x_t\right)d{x}_t{\eta }_{2,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_2\left(x_t\right)d{x}_t{\eta }_{2,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_2\left(x_t\right)d{x}_t{\eta }_{2,2}\left(0\right)\end{array}\right)+$$

$$+0.081·{\left[{\int }_{-\infty }^{+\infty }{\eta }_{2,1}^2\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\eta }_{2,2}^2\left(x\right)dx\right]}^{-1}\times$$

$$exp\left\{t{E}_0\right\}\left(\begin{array}{cc}{\int }_{-\infty }^{+\infty }{\phi }_0\left(x_t\right)d{x}_t{\eta }_{2,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\phi }_0\left(x_t\right)d{x}_t{\eta }_{2,2}\left(0\right)\\ {\int }_{-\infty }^{+\infty }{\psi }_0\left(x_t\right)d{x}_t{\eta }_{2,1}\left(0\right)& {\int }_{-\infty }^{+\infty }{\psi }_0\left(x_t\right)d{x}_t{\eta }_{2,2}\left(0\right)\end{array}\right).$$

For $s=0$, $t=1$, $a=1$, $b=1$, $A=5$, and $N=90$

$$E_2\approx -3.379,{\int }_{-\infty }^{+\infty }{\psi }_2\left(x_t\right)d{x}_t\approx -1.1035,{\int }_{-\infty }^{+\infty }{\phi }_2\left(x_t\right)d{x}_t\approx 0,$$

$${\psi }_2\left(0\right)\approx 0.6421,{\phi }_2\left(0\right)\approx -0.0849,$$

$${\int }_{-\infty }^{+\infty }{\eta }_{2,1}^2\left(x\right)dx+{\int }_{-\infty }^{+\infty }{\eta }_{2,2}^2\left(x\right)dx\approx 1.2655,{\eta }_{2,1}\left(0\right)\approx -0.0873,{\eta }_{2,2}\left(0\right)\approx 0.7108.$$

Thus, the term with $j=2$ contributes

$$\left(\begin{array}{cc}0& 0\\ 0.0007& -0.006\end{array}\right).$$

Taking into account the contribution of three terms, we obtain

$$\int exp\left\{{\int }_s^t\left[\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right]d\tau \right\}d\mu \left(x\right)\approx \left(\begin{array}{cc}0& 0\\ -0.0080& 0.2452\end{array}\right).$$

Statement 5 implies that for the elements of the matrix

$$exp\int exp\left\{{\int }_s^t\left[\beta \mathsf{\Phi }\left(x\left(\tau \right)\right)+V\left(x\left(\tau \right)\right)\right]d\tau \right\}d\mu \left(x\right)$$

with indices $(1,2)$ and $(2,1)$ for $V\left(x\right)=-x^2/2$, $\mathsf{\Phi }\left(x\right)=x^2/2$, the exact values are known and they are equal to zero. For comparison, Figure 29 shows the approximate values of the elements of the matrix valued integral with index $(2,1)$ for different values N.

Figure 29 shows that the approximate values of the elements of the matrix valued integral tend to the exact value as the number of intervals N increases into which the interval $[-A,A]$ is divided for the approximate calculation of the eigenvalues and eigenfunctions.

4. Conclusions

The method of functional integration finds ever wider application in various fields of science. Therefore, it becomes necessary to develop methods for the approximate calculation of functional integrals and generalize existing and newly obtained methods to new types of functional integrals. The methods we have proposed for calculating functional integrals generated by nonrelativistic Hamiltonians complement existing methods. As noted, formulas for a given degree of accuracy are especially effective for functional integrals over the space of functions defined on intervals of small length. For integrals determined over large time intervals, it is necessary to take formulas with a high degree of accuracy, which leads to high computational costs.

The efficiency of existing methods for the approximate calculation of functional integrals based on discretization of space and time also decreases with increasing time intervals. The method for calculating integrals based on expansion in terms of eigenfunctions of the Hamiltonian is effective for functional integrals defined over large time intervals. At present, the theory of approximate calculation of functional integrals generated by nonrelativistic Hamiltonians is well developed. The theory of approximate calculation of functional integrals generated by relativistic Hamiltonians is less developed. Therefore, generalizations of methods based on the approximation by polynomials and methods based on expansion in terms of eigenfunctions of a Hamiltonian to the case of relativistic Hamiltonians are very useful and desirable.

In further studies, it is planned to generalize these methods to functional integrals generated by stochastic models. In particular, it is planned to generalize the existing methods to functional integrals arising from the study of one-step processes by combinatorial and operator methods.