SOLUTION BY THE METHOD OF LEAST SQUARES OF A SPECTRAL PROBLEM FOR THE TWO-DIMENSIONAL HELMHOLTZ EQUATION WITH A BOUNDARY CONDITION

B. F. Krutoi

(Presented by the scientific seminar of the Department of Engineering and Computational Mathematics)

Suppose that in a two-dimensional, in general multiply connected, closed domain \(D^{(2)}=D^{(2)}+\Gamma\) the Helmholtz equation is given

\[ [\Delta + k^2]u(x,y)=\left[\frac{\partial^2}{\partial x^2}+\frac{\partial^2}{\partial y^2}+k^2\right]u(x,y)=f(x,y), \tag{1} \]

with a boundary condition of mixed type

\[ \left[p(x,y)+q(x,y)\frac{\partial}{\partial n}\right]u(x,y)\bigg|_{\Gamma}=h(x,y). \tag{2} \]

on the boundary \(\Gamma\) (Fig. 1 in [1]). Here \(k^2\) is a constant, \(\dfrac{\partial}{\partial n}u(x,y)\) is the derivative along the unit normal \(\bar n\) to the boundary \(\Gamma\), and \(f(x,y)\) and \(h(x,y)\) are given functions; moreover, since \(|\bar n|=1\), we have

\[ \begin{gathered} 1)\ \frac{\partial u}{\partial n} =\frac{\partial u}{\partial x}\frac{\partial x}{\partial n} +\frac{\partial u}{\partial y}\frac{\partial y}{\partial n}; \qquad 2)\ \frac{\partial x}{\partial n}=\cos(x,n)=\frac{\partial u}{\partial x}:N, \\[4pt] 3)\ \frac{\partial y}{\partial n}=\cos(y,n)=\frac{\partial u}{\partial y}:N; \qquad 4)\ N^2=\left(\frac{\partial u}{\partial x^2}\right)^2+ \left(\frac{\partial u}{\partial y}\right)^2 . \end{gathered} \tag{3} \]

Then, in view of the linearity of the operator \([\Delta=k^2]\), for the Helmholtz equation (1) with the boundary condition (2) there may be posed the corresponding spectral problem:

\[ \left\{ \begin{array}{l} \text{determine the sets (spectrum) of eigenvalues and eigenfunctions }\\ w(x,y)\text{ of the operator }[\Delta+k^2],\text{ satisfying in the domain}\\ \bar D^{(2)}=D^{(2)}+\Gamma\text{ the requirement} \end{array} \right. \]

\[ [\Delta+k^2]w(x,y)=\lambda w(x,y) \tag{4} \]

and the homogeneous boundary condition

\[ \left[p(x,y)+q(x,y)\frac{\partial}{\partial n}\right]w(x,y)\bigg|_{\Gamma}=0, \tag{5} \]

on the boundary \(\Gamma\). A solution of the indicated problem by the Fourier method set forth in [2] is impossible, since the coefficients \(p(x,y)\) and \(q(x,y)\) of the boundary condition (5) depend on both \(x\) and \(y\). Therefore we shall solve the given problem by the method of least squares, which also makes it possible to estimate the accuracy of the solution obtained at any point of the domain \(\overline{D}^{(2)}=D^{(2)}+\Gamma\). The essence of the proposed solution, which we shall divide into several successive stages, is as follows.

1. We represent the required function \(w(x,y)\) by its approximate expression \(\varphi(x,y)\)—a finite expansion in \(r\) linearly independent trial functions \(\varphi_\nu(x,y)\), so that

\[ w(x,y)\approx \varphi(x,y)=\sum_{\nu=1}^{r} c_\nu \varphi_\nu(x,y), \tag{6} \]

where \(c_\nu\) are unknown coefficients.

Let us now prescribe on the boundary \(\Gamma\) some number \(\gamma\) of uniformly distributed points \(M_\Gamma=M(x_\Gamma,y_\Gamma)\), and in the domain \(D^{(2)}\) prescribe uniformly some number \(\sigma\) of points \(M_s=M(x_s,y_s)\) (Fig. 1 in [1]). Then, for the \(\gamma\) points \(M_\Gamma\), we compute the boundary condition (5), and for the \(\sigma\) points \(M_s\), we compute equality (4), proceeding from the approximate representation (6) for the required function \(w(x,y)\) and taking \(\gamma+\sigma=\omega>r\). As a result we obtain a system of \(\omega>r\) equations

\[ \left\{ \begin{aligned} 1)\;& \sum_{\nu=1}^{r} \left[ p(x_\Gamma,y_\Gamma)+q(x_\Gamma,y_\Gamma)\frac{\partial}{\partial n} \right]\varphi_\nu(x_\Gamma,y_\Gamma)c_\nu=\varepsilon_\Gamma, \\ &\hspace{4.5em}(\Gamma=1,2,\ldots,\gamma) \\[0.6em] 2)\;& \sum_{\nu=1}^{r} \left[ \left(\frac{\partial^2}{\partial x^2}+\frac{\partial^2}{\partial y^2}+k^2\right) \varphi_\nu(x_s,y_s)-\lambda\varphi_\nu(x_s,y_s) \right]c_\nu=\varepsilon_s, \\ &\hspace{4.5em}(s=1,2,\ldots,\sigma) \end{aligned} \right. \tag{7} \]

in which \(\varepsilon_\Gamma,\varepsilon_s\) are certain errors caused by the inaccuracy of representation (6) for the required function \(w(x,y)\).

Introduce the notation

\[ 1)\quad \left[ p(x_\Gamma,y_\Gamma)+q(x_\Gamma,y_\Gamma)\frac{\partial}{\partial n} \right]\varphi_\nu(x_\Gamma,y_\Gamma)=a_{\Gamma\nu}=a_{\nu\Gamma}^{*}, \]

\[ 2)\quad \left(\frac{\partial^2}{\partial x^2}+\frac{\partial^2}{\partial y^2}+k^2\right) \varphi_\nu(x_s,y_s)=b_{s\nu}=b_{\nu s}^{*}, \tag{8} \]

\[ 3)\quad \varphi_\nu(x_s,y_s)=\varphi_{s\nu}=\varphi_{\nu s}^{*}. \]

Then we reduce the system (7) of \(\gamma+\sigma=\omega\) equations to the more compact form

\[ \left\{ \begin{aligned} 1)\;& \sum_{\nu=1}^{r} a_{\Gamma\nu}c_\nu=\varepsilon_\Gamma, \qquad (\Gamma=1,2,\ldots,\gamma) \\[0.6em] 2)\;& \sum_{\nu=1}^{r}(b_{s\nu}-\lambda\varphi_{s\nu})c_\nu=\varepsilon_s, \qquad (s=1,2,\ldots,\sigma). \end{aligned} \right. \tag{9} \]

The partial system (9.1) consists of \(\gamma\) algebraic equations of the first degree. The partial system (9.2) consists of \(\sigma\) algebraic equations

of the second degree, since it contains terms \(\lambda \varphi_{s\nu}c_\nu\), where \(\lambda\) and \(c_\nu\) are unknown. Therefore, in view of the solution of system (9) by least squares from the equations \(\nu+\sigma=\omega>r\), we must first find approximate values \(\lambda^0, c_\nu^0\) for the indicated unknowns \(\lambda, c_\nu\). At the same time one can compute the corresponding approximate values \(w^0(x,y)\) for the eigenfunctions \(w(x,y)\) of the operator \([\Delta+k^2]\).

1. Assuming that \(\sigma>r\), we shall compute \(\lambda^0, c_\nu^0\) as follows. Select any \(r\) equations from the partial system (9.2), so that in the end we obtain the system

\[ \sum_{\nu=1}^{r}(b_{\varkappa\nu}-\lambda\varphi_{\varkappa\nu})c_\nu=0, \qquad (\varkappa=1,2,\ldots,r), \tag{10} \]

of \(r\) algebraic equations of the second degree without free terms and with \(r+1\) unknowns \(\lambda, c_1, c_2,\ldots,c_r\). Bearing in mind the approximate representation (6) for \(w(x,y)\), we set ourselves the task of finding from system (10) such values of the \(r\) coefficients \(c_\nu\) in (6) that

\[ |C_1|+|C_2|+\cdots+|C_r|\ne 0. \tag{11} \]

Since the system of \(r\) equations (10) is homogeneous with respect to the \(r\) unknowns \(c_\nu\), requirement (11) will be fulfilled if the determinant \(b_{\varkappa\nu}-\lambda\varphi_{\varkappa\nu}\big|_1^r\) of system (10) is equal to zero,

\[ \big|b_{\varkappa\nu}-\lambda\varphi_{\varkappa\nu}\big|_1^r=0, \qquad (\varkappa,\nu=1,2,\ldots,r). \tag{12} \]

But the determinant \(\big|b_{\varkappa\nu}-\lambda\varphi_{\varkappa\nu}\big|_1^r\), after expansion, is a polynomial of degree \(r\) with respect to \(\lambda\). Consequently, condition (12) can also be written as

\[ \big|b_{\varkappa\nu}-\lambda\varphi_{\varkappa\nu}\big|_1^r = \sum_{s=0}^{r} g_s\lambda^s=0, \tag{13} \]

that is, the condition gives an algebraic equation of degree \(r\) with respect to the unknown \(\lambda\).

Solving equation (13), we find its \(r\) roots \(\lambda_s\), some of which may be multiple. Therefore the number \(n\) of distinct roots \(\lambda_\beta\) may be less than \(r\): \(n\le r\). Since these distinct values \(\lambda_\beta\) of the unknown \(\lambda\) have been obtained from the partial system of equations (10), we shall regard them as approximate and therefore denote them by \(\lambda_\beta^0\).

Next substitute into the system of equations (10), which contains no free terms, the \(n\) distinct roots \(\lambda_\beta\) of equation (13) that have been found. Denoting by \(c_{\nu\beta}^0=c_{1\beta}^0,c_{2\beta}^0,\ldots,c_{r\beta}^0\) the values of the unknowns \(c_\nu\) corresponding to the given \(\lambda_{0\beta}\), we write the result of substituting the \(n\) distinct roots \(\lambda_\beta\) into system (10) as follows:

\[ \text{a) } \sum_{\nu=1}^{r}(b_{\varkappa\nu}-\lambda_\beta^0\varphi_{\varkappa\nu})c_{\nu\beta}^0=0 \quad \text{or} \quad \text{b) } \sum_{\nu=1}^{r} b_{\varkappa\nu}c_{\nu\beta}^0 = \lambda_\beta^0\sum_{\nu=1}^{r}\varphi_{\varkappa\nu}c_{\nu\beta}^0, \tag{14} \]

\[ (\varkappa=1,2,\ldots,r;\qquad \beta=1,2,\ldots,n\le r), \]

which we shall also represent in matrix form as

\[ \text{a) } (b_{(rr)}-\lambda_\beta^0\varphi_{(rr)})c_{(r1)\cdot\beta}^0=0_{(r1)} \quad \text{or} \quad \text{b) } b_{(rr)}c_{(r1)\cdot\beta}^0 = \lambda_\beta^0\varphi_{(rr)}c_{(r1)\cdot\beta}^0. \tag{14*} \]

\[ (\beta=1,2,\ldots,n\le r). \]

It follows from this that, as a result of the indicated substitution, we have obtained a bundle of \(n \le r\) systems of equations of the form (10), and in each of these systems, instead of the unknown \(\lambda\), one of its different approximate values \(\lambda^0_\beta \in \lambda^0_1,\lambda^0_2,\ldots,\lambda^0_n\) has been substituted.

Since each of the \(n\) systems of equations entering into the bundle (14) consists of \(r\) homogeneous equations with respect to the corresponding \(r\) unknowns \(c^0_{\nu\beta}=c^0_{1\beta},c^0_{2\beta},\ldots,c^0_{r\beta}\), one of these unknowns remains arbitrary. Therefore we shall assume that in each system of homogeneous equations from the bundle (14), for the given \(\lambda^0_\beta\), the last unknown \(c^0_{r\beta}=1\). Then the remaining \(r-1\) unknowns \(c^0_{1\beta},c^0_{2\beta},\ldots,c^0_{r-1,\beta}\) for the given \(\lambda^0_\beta\) are determined completely uniquely by their \(\beta\)-th system of \(r-1\) equations, so that instead of (14) we shall have

\[ \sum_{\nu=1}^{r-1}(b_{\chi\nu}-\lambda^0_\beta\varphi_{\chi\nu})c^0_{\nu\beta} +(b_{\chi r}-\lambda^0_\beta\varphi_{\chi r})=0,\quad c^0_{r\beta}=1; \tag{15} \]

\[ (\chi=1,2,\ldots,r-1;\ \beta=1,2,\ldots,n\le r). \]

The set of \(r\times n\) numbers \(c^0_{\nu\beta}\) obtained by the method described above forms, in general, a rectangular matrix \(c^0_{(rn)}\).

\[ c^0_{(rn)}= \left[ \begin{array}{cccccc} c^0_{11} & c^0_{12} & \cdots & c^0_{1\beta} & \cdots & c^0_{1n}\\ c^0_{21} & c^0_{22} & \cdots & c^0_{2\beta} & \cdots & c^0_{2n}\\ \cdot & \cdot & \cdots & \cdot & \cdots & \cdot\\ c^0_{\nu1} & c^0_{\nu2} & \cdots & c^0_{\nu\beta} & \cdots & c^0_{\nu n}\\ \cdot & \cdot & \cdots & \cdot & \cdots & \cdot\\ c^0_{r1} & c^0_{r2} & \cdots & c^0_{r\beta} & \cdots & c^0_{rn} \end{array} \right],\quad c^0_{r\beta}=1, \tag{16} \]

which becomes a square matrix \(c^0_{(rr)}\) when \(n=r\). If each column \(c^0_{1\beta}, c^0_{2\beta},\ldots,c^0_{\nu\beta},\ldots,c^0_{r\beta}\) of the matrix \(c^0_{(rn)}\) is regarded as a vector \(\bar c^0_\beta\), then this matrix is written more compactly as follows:

\[ 1)\quad \bar c^0_\beta= \left[ \begin{array}{c} c^0_{1\beta}\\ c^0_{2\beta}\\ \cdot\\ \cdot\\ \cdot\\ c^0_{\nu\beta}\\ \cdot\\ \cdot\\ \cdot\\ c^0_{r\beta} \end{array} \right] =c^0_{(r1)\cdot\beta}, \quad 2)\quad c^0_{(rn)}=(\bar c^0_1,\bar c^0_2,\ldots,\bar c^0_\beta,\ldots,\bar c^0_n). \tag{16a} \]

On this basis the bundle (14) of \(n\le r\) systems of equations (10) can also be represented in the following vector form:

\[ \text{a) }\ (B-\lambda_\beta\Phi)\bar c^0_\beta=\bar 0_\beta,\qquad \text{b) }\ \bar B\bar c^0_\beta=\lambda_\beta\Phi\bar c^0_\beta, \tag{14**} \]

\[ (\beta=1,2,\ldots,n\le r) \]

where \(B\) and \(\Phi\) are linear operators corresponding to the matrices \(b_{(rr)}\), \(\varphi_{(rr)}\) in some reference basis,

1. Suppose further that from the solution of equation (13) we have obtained exactly \(r\) approximate eigenvalues \(\lambda_\chi^0 \in \lambda_1^0,\lambda_2^0,\ldots,\lambda_r^0\) of the operator \([\Delta+k^2]\). Substituting then the \(r\) different \(\lambda_\chi^0\) found into relation (14), we determine from it \(r\) vectors \(\bar c_\chi^0 \in \bar c_1^0,\bar c_2^0,\ldots,\bar c_r^0\) belonging to these \(\lambda_\chi^0\). Then equality (6) will give us exactly \(r\) approximate representations \(w_\chi^0(x,y)\) for the eigenfunctions \(w(x,y)\) of the operator \([\Delta+k^2]\), so that instead of (6) we obtain the following refined expression:

\[ w(x,y)\in w_\chi(x,y)\simeq w_\chi^0(x,y)= \sum_{\nu=1}^{r}\dot c_{\chi,\nu}^{0}\psi_\nu(x,y), \tag{17} \]

\[ (\chi=1,2,\ldots,r). \]

1. Suppose now that from equation (13) we have obtained only \(n<r\) different approximate eigenvalues \(\lambda_\beta^0 \in \lambda_1^0,\lambda_2^0,\ldots,\lambda_n^0\) of the operator \((\Delta+k^2)\). If we substitute the \(n\) different \(\lambda_\beta^0\) found into relation (14), then from it we determine only \(n\) corresponding vectors \(\bar c_\beta^0=\bar c_1^0,\bar c_2^0,\ldots,\bar c_n^0\). For substitution into representation (6) for the eigenfunctions \(w(x,y)\), however, we need exactly \(r\) independent vectors \(\bar c_\chi^0\). Consequently, another \(r-n\) independent vectors \(\bar g_\beta^0\) are lacking; we shall call them adjoined and form them in the following way.

First of all, it is clear that for each \(\partial_\beta\)-fold root \(\lambda_\beta\) of equation (13) we can construct, by the method described above, only one vector \(\bar c_\beta^0\) instead of the \(\partial_\beta\) necessary independent vectors, which we denote by \(\bar c_{\beta\alpha}^0\in \bar c_{\beta1}^0,\bar c_{\beta2}^0,\ldots,\bar c_{\beta\cdot\partial_\beta}^0\). If we set \(\bar c_{\beta1}^0=\bar c_\beta^0\), then all \(\partial_\beta\) independent vectors \(c_{\beta\alpha}^0\) for the given \(\lambda_\beta\) can be constructed in many ways. For example, as the solution of a system of \(\partial_\beta\) independent vector equations of the following special form:

\[ \left\{ \begin{aligned} B\bar c_{\beta1}^{0}&=\lambda_\beta^0\Phi \bar c_{\beta1}^{0} +0\bar c_{\beta2}^{0}+0\bar c_{\beta3}^{0}+\cdots +0\bar c_{\beta\cdot\partial_\beta-1}^{0} +0\bar c_{\beta\cdot\partial_\beta}^{0},\\ B\bar c_{\beta2}^{0}&=1\bar c_{\beta1}^{0} +\lambda_\beta^0\Phi \bar c_{\beta2}^{0} +0\bar c_{\beta3}^{0}+\cdots +0\bar c_{\beta\cdot\partial_\beta-1}^{0} +0\bar c_{\beta\cdot\partial_\beta}^{0},\\ B\bar c_{\beta3}^{0}&=0\bar c_{\beta1}^{0} +1\bar c_{\beta2}^{0} +\lambda_\beta^0\Phi \bar c_{\beta3}^{0} +\cdots +0\bar c_{\beta\cdot\partial_\beta-1}^{0} +0\bar c_{\beta\cdot\partial_\beta}^{0},\\ &\cdot\quad \cdot\quad \cdot\quad \cdot\quad \cdot\quad \cdot\quad \cdot\quad \cdot\\ B\bar c_{\beta\cdot\partial_\beta}^{0}&=0\bar c_{\beta1}^{0} +0\bar c_{\beta2}^{0} +0\bar c_{\beta3}^{0}+\cdots +1\bar c_{\beta\cdot\partial_\beta-1}^{0} +\lambda_\beta^0\Phi \bar c_{\beta\cdot\partial_\beta}^{0}, \end{aligned} \right. \tag{18} \]

where \(B\) and \(\Phi\) are linear operators corresponding to the matrices \(b_{(rr)}\) and \(\varphi_{(rr)}\).

The solution of the first of the vector equations (19), which is rather difficult when the number of components \(c_{\beta1,\nu}^{0}\), \(\nu=1,2,\ldots,r\), is large, can be carried out by a generalized rotation method and will be presented in a separate paper. But if the vector \(\bar c_{\beta1}^{0}\) is found, then the remaining \(\partial_\beta-1\) vectors \(\bar c_{\beta\alpha}^{0}\) are determined, obviously, as follows:

\[ \bar c_{\beta\alpha}^{0}=(B-\lambda_\beta^0\Phi)^{\alpha-1}\bar c_{\beta1}^{0} =(B-\lambda_\beta\Phi)^{\alpha-1}\bar c_\beta^{0}. \tag{19} \]

\[ (\alpha=2,3,\ldots,\partial_\beta). \]

Having carried out, according to (18) and (19), the determination of all \(r=\partial_1+\partial_2+\cdots,\partial_\beta+\cdots,\partial_n\) principal \(c_\beta^0=c_{\beta1}^0\) and adjoined \(c_{\beta\alpha}^0=c_{\beta2}^0,c_{\beta3}^0,\ldots,c_{\beta\cdot\partial_\beta}^0\)

vectors and substituting them into the representation (6) for the eigenfunctions \(w'(x,y)\) of the operator \([\Delta+k^2]\), we obtain a computational expression for a finite set \(w_x^0(x,y)\approx w(x,y)\) of the same form (17) as before.

However, if expression (17) is written in matrix form

\[ w_{(r1)}(x,y)=\overset{*0}{c}_{(rr)}\varphi_{(r1)}(x,y) \tag{17*} \]

and it is taken into account that

\[ c^0_{(rr)}=\left[\,c^0_{(rn)}\ c^0_{(r.r-n)}\,\right] = \left[\,c^0_{(rn)}\ f^0_{(r.r-n)}\,\right], \]

where \(c^0_{(rn)}\) and \(c^0_{(r.r-n)}=f^0_{(r.r-n)}\) are matrices of the principal \(c^0_{\nu\beta}\) and adjoined \(c^0_{\nu\gamma}=f^0_{\nu\gamma}\), \((\gamma=n+1,n+2,\ldots,r)\) coefficients, then equality (17) is written in more expanded form as follows:

\[ w_{(r1)}(x,y)= \begin{bmatrix} \overset{*0}{c}_{(nr)}\\ \overset{*0}{f}_{(r-n.r)} \end{bmatrix} \varphi_{(r1)}(x,y) = \begin{bmatrix} c^0_{(nr)}\varphi_{(r1)}\\ \overset{*0}{f}_{(r-n.r)}\varphi_{(r1)} \end{bmatrix}. \tag{17*.1} \]

Hence the following refined representation for \(w_x(x,y)\) follows:

\[ \begin{cases} w_\beta(x,y)=\displaystyle\sum_{\nu=1}^{r}\overset{*0}{c}_{\beta\nu}\varphi_\nu(x,y),\\[6pt] \qquad (x=1,2,\ldots,\beta,\ldots,n)\\[6pt] w_\gamma(x,y)=\displaystyle\sum_{\nu=1}^{r}\overset{*0}{f}_{\gamma\nu}\varphi_\nu(x,y)\\[6pt] \qquad (x=n+1,n+2,\ldots,\gamma,\ldots,r). \end{cases} \tag{17.1} \]

1. Let us now rewrite the system of error equations (9) in a more expanded form, following from the refined representations

\[ 1)\ \lambda\in\lambda_1,\lambda_2,\ldots,\lambda_\beta,\ldots,\lambda_n,\quad 2)\ c_\nu=c_{\nu x},\quad (\nu,x=1,2,\ldots,r) \tag{20} \]

for the unknowns \(\lambda\) and \(c_\nu\). Then the system (9) is written more fully as follows:

\[ \begin{cases} 1)\ \displaystyle\sum_{\nu=1}^{r} a_{\Gamma\nu}c_{\nu\beta}=\varepsilon_{\Gamma\beta}\\[10pt] \qquad \begin{pmatrix} \Gamma=1,2,\ldots,\gamma\\ s=1,2,\ldots,\sigma \end{pmatrix}\\[10pt] 2)\ \displaystyle\sum_{\nu=1}^{r}\left(b_{s\nu}c_{\nu\beta}-\varphi_{s\nu}c_{\nu\beta}\lambda_\beta\right)=\hat{\varepsilon}_{s\beta},\\[8pt] (\beta=1,2,\ldots,n\le r;\ \sigma>r,\ \gamma+\sigma=\omega>r+1) \end{cases} \tag{21} \]

or, in matrix form,

\[ \begin{cases} 1)\ a_{(\gamma r)}c_{(rn)}=\varepsilon_{(\gamma n)}\\ 2)\ b_{(\sigma r)}c_{(rn)}-\varphi_{(\sigma r)}c_{(rn)}\lambda_{(nn)}=\varepsilon_{(\sigma n)}. \end{cases} \qquad (n\le r;\ \sigma>n,\ \gamma+\sigma=\omega>r+1). \tag{21*} \]

Thus, the system of error equations (9) is in fact the bundle (21) of \(2n\) systems of equations; moreover, for a given \(\lambda_\beta\), the system (21.1) contains \(\gamma\) equations, and the system (21.2) contains \(\sigma>r\) equations. Consequently, \(\gamma+\sigma>r+1\), that is, for each \(\lambda_\beta\) the number of equations (21), equal to \(\gamma+\sigma=\omega\), will be greater than the number \(m+1\) co-

corresponding unknowns \(\lambda_\beta, c_{\nu\beta}\). Therefore the solution of the system (21*) consisting of \(2n\) sets of equations must be carried out by the method of least squares, which also makes it possible to estimate reliably the accuracy of the eigenvalues \(\lambda_\beta\) determined from these equations and of the eigenvectors \(\overline{c}_\beta=\overline{c}_{\beta 1}\) of the operator \([\Delta+k^2]\). In addition, we shall then be able to estimate the accuracy of the associated vectors \(c_{\nu\alpha}\), \(\alpha=2,3,\ldots,\delta_\beta\), determined according to (19), for eigenvalues \(\lambda_\beta\) of multiplicity \(\delta_\beta\). In conclusion we shall estimate, by the method of least squares, the accuracy of the finite set \(w_\alpha(x,y)\), determined according to (17), of eigenfunctions \(w(x,y)\) of the operator \([\Delta+k^2]\).

1. Proceeding to solve by the method of least squares the system (21) consisting of \(2n\) sets of equations, we represent the unknowns \(\lambda_\beta, c_{\nu\beta}\) in (21) in the following form:

\[ 1)\ \lambda_\beta=\lambda_\beta^{0}+\delta\lambda_\beta,\qquad 2)\ c_{\nu\beta}=c_{\nu\beta}^{0}+\delta c_{\nu\beta}, \tag{22} \]

where \(\lambda_\beta^{0}, c_{\nu\beta}^{0}\) are previously found approximate, but sufficiently accurate, values for \(\lambda_\beta, c_{\nu\beta}\). Therefore the corrections \(\delta\lambda_\beta, \delta c_{\nu\beta}\) will be regarded as so small that their squares \((\delta\lambda_\beta)^2\), \((\delta c_{\nu\beta})^2\) and the products \(\delta\lambda_\beta, \delta c_{\nu\beta}\) may be neglected. Substituting then the representations (22) for \(\lambda_\beta, c_{\nu\beta}\) into equations (21), and taking into account only the remarks just made, we obtain

\[ \left\{ \begin{aligned} 1)\ &\sum_{\nu=1}^{r} a_{r\nu}\,\delta c_{\nu\beta}+\zh_{r\beta}=\varepsilon_{r\beta},\\ 2)\ &\sum_{\nu=1}^{r}\left[b_{s\nu}\,\delta c_{\nu\beta} -(\varphi_{s\nu}c_{\nu\beta}^{0})\,\delta\lambda_\beta -(\lambda_\beta^{0}\varphi_{s\nu})\,\delta c_{\nu\beta}\right] +w_{s\beta}=\varepsilon_{s\beta}, \end{aligned} \right. \tag{23} \]

\[ (r=1,2,\ldots,\gamma;\ s=1,2,\ldots,\sigma) \]

\[ (\beta=1,2,\ldots,n\le r;\ \sigma>r,\ \gamma+\sigma>m+1) \]

where the notation is

\[ 1)\ \sum_{\nu=1}^{r}a_{r\nu}c_{\nu\beta}^{0}=\zh_{r\beta},\qquad 2)\ \sum_{\nu=1}^{r}\left[b_{s\nu}c_{\nu\beta}^{0} -(\varphi_{s\nu}c_{\nu\beta}^{0})\lambda_\beta^{0}\right]=-w_{s\beta}. \tag{24} \]

Introduce the additional notation

\[ 1)\ \delta\lambda_\beta=\eta_\beta,\qquad 2)\ \delta c_{\nu\beta}=\xi_{\nu\beta},\qquad 3)\ \sum_{\nu=1}^{r}\varphi_{s\nu}c_{\nu\beta}^{0}=-p_{s\beta}, \tag{25} \]

\[ 4)\ \lambda_\beta^{0}\varphi_{s\nu}=q_{s\nu}^{(\beta)},\qquad 5)\ b_{s\nu}-q_{s\nu}^{(\beta)}=g_{s\nu}^{(\beta)}. \]

Then the collection of equations (23) will be written in the following compressed form:

\[ \left\{ \begin{aligned} 1)\ &\sum_{\nu=1}^{r}a_{r\nu}\xi_{\nu\beta}+\zh_{r\beta}=\varepsilon_{r\beta},\\ 2)\ &\sum_{\nu=1}^{r}g_{s\nu}^{(\beta)}\xi_{\nu\beta} +p_{s\beta}\eta_\beta+w_{s\beta}=\varepsilon_{s\beta}, \end{aligned} \right. \tag{26} \]

\[ (r=1,2,\ldots,\gamma;\ s=1,2,\ldots,\sigma) \]

\[ (\beta=1,2,\ldots,n\le r;\ \sigma<r,\ \gamma+\sigma=\omega>r+1) \]

or, in matrix form,

\[ \begin{aligned} 1)\quad & a_{(\gamma r)}\xi_{(rn)}+\mathcal{H}_{(\gamma n)}=\varepsilon_{(\gamma n)},\\ 2)\quad & g_{(\sigma r)}\xi_{(rn)}+p_{(\sigma n)}\eta_{(nn)}+\mathcal{W}_{(\sigma n)}=\varepsilon_{(\sigma n)}, \end{aligned} \tag{26*} \]

\[ (n\le r;\ \sigma<r,\ \gamma+\sigma=\omega>r+1), \]

where \(\eta_{(nn)}\) is a diagonal matrix. If we now introduce in (26*) the notations

\[ \begin{aligned} 1)\quad& \begin{bmatrix} a_{(\gamma r)}\\ g_{(\sigma r)} \end{bmatrix} =A_{(\gamma+\sigma,r)}=A_{(\omega r)},\qquad 2)\quad& \begin{bmatrix} 0_{(\gamma n)}\\ p_{(\sigma n)} \end{bmatrix} =B_{(\gamma+\sigma,n)}=B_{(\omega n)},\\ 3)\quad& \begin{bmatrix} \mathcal{H}_{(\gamma n)}\\ \mathcal{W}_{(\sigma n)} \end{bmatrix} =K_{(\gamma+\sigma,n)}=K_{(\omega n)},\qquad 4)\quad& \begin{bmatrix} \varepsilon_{(\gamma n)}\\ \varepsilon_{(\sigma n)} \end{bmatrix} =\varepsilon_{(\gamma+\sigma,n)}=\varepsilon_{(\omega n)}, \end{aligned} \tag{27} \]

then the bundle (21) of \(2n\) systems of equations takes the following final form:

\[ A_{(\omega r)}\xi_{(rn)}+B_{(\omega n)}\eta_{(nn)}+K_{(\omega n)}=\varepsilon_{(\omega n)},\quad(\omega>n). \tag{21—28} \]

Since in the matrix equation (21—28) it is assumed that \(\omega\cdot r>r\cdot n\) and \(\omega\cdot n>n\cdot n\), hence \(\omega>n\), we solve this matrix equation under the Gauss condition

\[ \sum_{\beta=1}^{r}\varepsilon^{(\beta)}_{(1\omega)}\varepsilon^{(\beta)}_{(\omega1)}=Q^{2}=\text{min.}, \tag{29} \]

whence two matrix equations follow [3]:

\[ \begin{aligned} &A^{*}_{(r\omega)}A_{(\omega r)}\xi_{(rn)} +A^{*}_{(r\omega)}B_{(\omega n)}\eta_{(nn)} +A^{*}_{(r\omega)}K_{(\omega n)}=0_{(rn)},\\ &B^{*}_{(n\omega)}A_{(\omega r)}\xi_{(rn)} +B^{*}_{(n\omega)}B_{(\omega n)}\eta_{(nn)} +B^{*}_{(n\omega)}K_{(\omega n)}=0_{(nn)} \end{aligned} \tag{30} \]

for determining the two unknown matrices \(\xi_{(rn)}\) and \(\eta_{(nn)}\). If we introduce the notations

\[ \begin{aligned} 1)\quad& A^{*}_{(r\omega)}A_{(\omega r)}=P_{(rr)},\qquad 2)\quad A^{*}_{(r\omega)}B_{(\omega n)}=\bigl(B^{*}_{(n\omega)}A_{(\omega r)}\bigr)^{*}=D_{(nr)},\\ 3)\quad& B^{*}_{(n\omega)}B_{(\omega n)}=R_{(nn)},\qquad 4)\quad A^{*}_{(r\omega)}K_{(\omega n)}=S_{(rn)},\qquad 5)\quad B^{*}_{(n\omega)}K_{(\omega n)}=T_{(nn)}, \end{aligned} \tag{31} \]

then the matrix equations (30) are written more compactly as

\[ \begin{cases} P_{(rr)}\xi_{(rn)}+D_{(rn)}\eta_{(nn)}+S_{(rn)}=0_{(rn)},\\ D^{*}_{(nr)}\xi_{(rn)}+R_{(nn)}\eta_{(nn)}+T_{(nn)}=0_{(nn)}. \end{cases} \tag{32} \]

Having found from this the correction matrices \(\xi_{(rn)}\) and \(\eta_{(nn)}\), we then compute, according to (25) and (22), the required principal matrices \(\lambda_{(nn)}\) and \(c_{(rn)}\):

\[ 1)\quad \lambda_{(nn)}=\lambda^{0}_{(nn)}+\eta_{(nn)} =\lambda^{0}_{(nn)}+\delta\lambda_{(nn)},\qquad 2)\quad c_{(rn)}=c^{0}_{(rn)}+\xi_{(rn)} =c^{0}_{(rn)}+\delta c_{(rn)}. \]

This completes the determination, by the least-squares method, of the principal matrices \(\lambda_{(nn)}\) and \(c_{(rn)}\).

1. Having computed the matrix \(c_{(rn)}\) of eigenvectors \(\overline c_{\beta}=c^{(\beta)}_{(r1)}\) for the bundle of equations (21*), we then compute, according to (19), the corresponding associated vectors \(\overline c_{\beta\alpha}\), \(\alpha=2,3,\ldots,\partial_{\beta}\), for those eigenvalues \(\lambda_{\beta}\) of the operator \([\Delta+k^{2}]\) which have multiplicity \(\partial_{\beta}\). Now we can write the already refined, by the least-squares method, ...

... a complete set of \(\omega_\chi(x,y)\), \((k=1,2,\ldots,r)\), of values for the eigenfunctions \(w(x,y)\) of the operator \((\Delta+k^2)\), which is obtained when, in (17), the approximate \(c_{(\chi\nu)}^{*0}\) are replaced by the refined \(c_{\chi\nu}^{*}\),

\[ w(x,y)\leftarrow w_\chi(x,y)=\sum_{\nu=1}^{r} c_{\chi\nu}^{*}\varphi_\nu(x,y). \tag{17a} \]

It is clear that in the case when all \(r\) roots \(\lambda_\beta\) of the equation \(|b_{(rr)}-\lambda\varphi_{(rr)}|=0\) are simple, the calculation according to (19) of all \(r-n\) adjoint vectors \(c_{\beta\alpha}\) is unnecessary, and immediately after solving the bundle of equations (21) we proceed to the calculation of the set of eigenfunctions \(w(x,y)\) according to (17a).

1. In conclusion, let us determine the mean-square spreads \(m^2(\lambda_\beta)\) and \(m^2[w_\chi(x,y)]\) of the eigenvalues \(\lambda_\beta\) and eigenfunctions \(w_\chi(x,y)\) of the operator \([\Delta+k^2]\) found by the least-squares method.

Generalizing the known rules for estimating accuracy in the least-squares method ([1], [3]), we first compute the mean-square spread \(\mu^2\) for the set of \(n\le r\) quantities \(\lambda_\mu\) and \(r\cdot n\le r\cdot r\) quantities \(c_{\nu\beta}\), specified by means of \(\gamma+\sigma=\omega\) error equations (21). We find that

\[ \mu^2= \frac{Q^2}{\omega n-(r n+n)} = \frac{\sum_{i=1}^{n}\varepsilon_{(1\omega)}^{(i)}\varepsilon_{(\omega 1)}^{(i)}} {\omega n-(r+1)n} = \frac{\sum_{\beta=1}^{r}\left(\sum_{\rho=1}^{\gamma}\varepsilon_{\rho\beta}^{2} +\sum_{s=1}^{\sigma}\varepsilon_{s\beta}^{2}\right)} {[\omega-(r+1)]n}. \tag{33} \]

Next we compute the mean-square spreads \(m^2(\lambda_\beta)\) and \(m^2(c_{\nu\beta})\), \((\nu=1,2,\ldots,r;\ \beta=1,2,\ldots,n\le r)\). For this purpose, let us rewrite the pair of matrix equations (32) in the form of a bundle of \(n\) pairs of equations corresponding to distinct \(\beta=1,2,\ldots,n\le r\):

\[ \begin{bmatrix} P_{(rr)}^{(\beta)} & D_{(rn)}^{(\beta)}\\ \ddot D_{(nr)}^{(\beta)} & R_{(nn)}^{(\beta)} \end{bmatrix} \begin{bmatrix} \xi_{(r1)}^{(\beta)}\\ \eta_{(n1)}^{(\beta)} \end{bmatrix} + \begin{bmatrix} S_{(r1)}^{(\beta)}\\ T_{(n1)}^{(\beta)} \end{bmatrix} = \begin{bmatrix} 0_{(r1)}^{(\beta)}\\ 0_{(n1)}^{(\beta)} \end{bmatrix}, \tag{32} \]

whence it follows that

\[ \begin{bmatrix} \xi_{(r1)}^{(\beta)}\\ \eta_{(n1)}^{(\beta)} \end{bmatrix} = - \begin{bmatrix} P_{(rr)}^{(\beta)} & D_{(rn)}^{(\beta)}\\ \ddot D_{(nr)}^{(\beta)} & R_{(nn)}^{(\beta)} \end{bmatrix}^{-1} \begin{bmatrix} S_{(r1)}^{(\beta)}\\ T_{(n1)}^{(\beta)} \end{bmatrix} = - \begin{bmatrix} \pi_{(rr)}^{(\beta)} & \hat\delta_{(rn)}^{(\beta)}\\ \ddot\delta_{(rn)}^{(\beta)} & \rho_{(nn)}^{(\beta)} \end{bmatrix} \begin{bmatrix} S_{(r1)}^{(\beta)}\\ T_{(n1)}^{(\beta)} \end{bmatrix}. \tag{33} \]

Then

\[ \text{1) }\ m^2(c_{\nu\beta})=\mu^2\pi_{\nu\nu}^{(\beta)},\qquad \text{2) }\ m^2(\lambda_i)=\mu^2\rho_{i=\beta\cdot r+\beta}^{(\beta)} . \tag{34} \]

\[ (\nu=1,2,\ldots,r;\ \beta=1,2,\ldots,n\le r) \]

Now we must calculate the mean-square spreads \(m^2(c_{\beta\alpha\delta})\), \((\delta=1,2,\ldots,r)\), of the components \(c_{\beta\alpha\delta}\), for the adjoint vectors \(\bar c_{\beta\alpha}\) belonging to some eigenvalue \(\lambda_\beta\) of multiplicity \(\partial_\beta\). According to (19),

\[ \bar c_{\beta\alpha}=(B-\lambda_\beta\Phi)^{\alpha-1}\bar c_{\beta1},\quad \bar c_{\beta1}=\bar c_\beta,\quad (\alpha=2,3,\ldots,\partial_\beta), \tag{19} \]

or, in matrix form,

\[ c_{(r1)}^{\beta\alpha} = \left(b_{(rr)}-\lambda_\beta\varphi_{(rr)}\right)^{\alpha-1} c_{(r1)}^{\beta1} = d_{(rr)}^{\beta\alpha}c_{(r1)}^{\beta}. \tag{19*} \]

and, in components,

\[ c_{\beta\alpha\cdot\delta}=c_{\delta}^{(\beta\alpha)} =\sum_{\nu=1}^{r} d_{\delta\nu}^{(\beta\alpha)}c_{\nu\beta},\quad (\delta=1,2,\ldots,r;\ \beta=1,2,\ldots,n\le r). \tag{19**} \]

Hence, by the general rule for calculating the mean-square scatter \(m_F^2\) of the function \(F(x_1,x_2,\ldots,x_j,\ldots,x_r)\) of \(r\) unknowns \(x_j\) connected by error equations, we obtain

\[ m^2(c_{\beta\alpha\cdot\delta})=\mu^2\sum_{\nu,\xi=1}^{r} (\pi_{\nu\xi}d_{\delta\nu}d_{\xi\delta})_{\beta\alpha}. \tag{35} \]

\[ \beta=1,2,\ldots,n\le r;\ \alpha=2,3,\ldots,\partial_\beta). \]

Finally, let us calculate the mean-square scatters \(m^2[w_\chi(x,y)]\) of the values \(w_\chi(x,y)\), found by the least-squares method, for the eigenfunctions \(w(x,y)\) of the operator \([\Delta+k^2]\). For this purpose, in the refined representation (17.1) for the eigenfunctions \(w_\chi(x,y)\) of the operator \([\Delta+k^2]\), we replace the approximate \(c_{\nu\beta}^{0}\) and \(c_{\beta\alpha\cdot\nu}^{0}=f_{\nu\gamma}^{0}\) by their adjusted values \(c_{\nu\beta}\) and \(f_{\nu\gamma}\). In addition, we take into account the expression (19) of the attached coefficients \(c_{\beta\alpha\cdot\delta}=f_{\delta\gamma}\) in terms of the principal coefficients \(c_{\nu\beta}\), which we write in a more convenient form

\[ \text{a) }\quad f_{\nu\gamma}=\sum_{\delta=1}^{r}d_{\nu\delta}c_{\delta\gamma},\qquad \text{b) }\quad f_{\gamma\nu}^{*}=\sum_{\delta=1}^{r}c_{\gamma\delta}^{*}d_{\delta\nu}^{*}, \tag{36} \]

\[ (\nu=1,2,\ldots,r;\ \gamma=n+1,n+2,\ldots,r). \]

Here \(c_{\delta\gamma}\) are the components of those principal vectors \(\overline{c_i}\), from their complete set \(\overline{c_\beta}\), \(\beta=1,2,\ldots,n\), which participated, according to (19), in the calculation of the corresponding attached vectors \(\overline{c_{\beta\alpha}}=\overline{f_\gamma}\).

We can now write the representation (17.1) for the eigenfunctions \(w_\chi(x,y)\) in the following expanded form:

\[ w_\chi(x,y)= \begin{cases} \displaystyle w_\beta(x,y)=\sum_{\nu=1}^{r}c_{\beta\nu}^{*}\varphi_\nu(x,y), & (\chi=1,2,\ldots,\beta,\ldots,n),\\[1.2em] \displaystyle w_\gamma(x,y)=\sum_{\delta=1}^{r}c_{\gamma\delta}^{*} \left(\sum_{\nu=1}^{r}d_{\delta\nu}^{*}\varphi_\nu(x,y)\right) =\sum_{\delta=1}^{r}c_{\gamma\delta}^{*}\psi_\delta(x,y), & (\chi=n+1,n+2,\ldots,\gamma,\ldots,r), \end{cases} \tag{17.2} \]

where, in calculating \(w_\gamma(x,y)\), the order of summation has been changed and the notation

\[ \psi_\delta(x,y)=\sum_{\nu=1}^{r}d_{\delta\nu}^{*}\varphi_\nu(x,y) \tag{37} \]

has been introduced.

From the obtained representation (17.2) for \(w_\chi(x,y)\) we find, by the general rule, that

\[ \text{a) }\quad m^2[w_\beta(x,y)]=\mu^2\varphi_{(1r)}^{*}(x,y)\pi_{(rr)}^{\beta}\varphi_{(r1)}(x,y), \]

\[ (\beta=1,2,\ldots,n\le r) \tag{38} \]

b)
\[ m^2[w_\gamma(x,y)] = \mu^2 \psi_{(rr)}^{*}(x,y)\pi_{(rr)}^\gamma \psi_{(r1)}(x,y). \]
\[ (\gamma=n+1,\ n+2,\ldots,r) \]

Here \(\pi_{(rr)}^\gamma\) is the value of the matrix \(\pi_{(rr)}\) for that eigenvalue \(\lambda_j \in \lambda_1,\lambda_2,\ldots,\lambda_n\) by means of which the associated vector \(C_{j\alpha}=\vec f_\gamma\) belonging to it was computed according to (19).

This completes the solution of the posed problem by the proposed method and the accuracy estimate for the solution obtained in this way. The method considered can, without any complications, be applied to solving and estimating the accuracy of a more general spectral problem

\[ 1)\quad L(x,y)u(x,y)=0,\qquad 2)\quad P(x,y)u(x,y)\big|_{\Gamma}=0 \tag{39} \]

with arbitrary linear differential operators \(L(x,y)\), \(P(x,y)\) and a function \(u(x,y)\) independent of time \(t\).

REFERENCES

1. B. F. Krutoi. Least-squares solution of a boundary-value problem for the two-dimensional Poisson equation. Izvestiya TPI, present collection.

2. V. S. Vladimirov. Equations of Mathematical Physics. Moscow, “Nauka,” 1971.

3. V. Iordan. Manual of Geodesy. Vol. 1. Editorial Bureau of GUGK, Moscow, 1939.