A Solution of the n-Dimensional Navier-Stokes System for an Incompressible Fluid with Cauchy Condition

1. Introduction

It is well known that, despite the nonlinearity of the incompressible NSE for a viscous fluid, numerous analytical studies have been conducted proposing various methods for obtaining solutions, typically under certain restrictions on the convective and viscous terms of the NSE [2,3,4,5,6,7,8], among others. For example, in specific cases, solutions have been obtained through the linearization of the equations of motion, neglecting the nonlinear convective terms under the assumption that they are sufficiently small compared to the frictional forces, particularly when the viscosity is relatively large (see, for example, Landau аnd Lifshitz [4], Schlichting and others). Under these assumptions, the resulting solutions to the NS system for an incompressible fluid show good agreement with experimental data, thereby confirming the general applicability of the NSE.

Therefore, from a physical standpoint, the derivation of the incompressible NSE is beyond the scope of this study, as there exists a vast body of fundamental work addressing these issues (see, for example, Caffarelli, Kohn and Nirenberg [2], Scheffer [3], Landau and Lifshith [4], Prantdl [5], Schlichting [6], among others). Moreover, in our case, NS problems involving sufficiently large viscosity are not considered, since such problems fall outside the scope of the Millennium Prize requirements (2000y., Fefferman [1]); on the other hand, to some extent, it can be said that the theory of such NS systems has already been developed.

Accordingly, in this paper, we do not attempt to review the extensive literature on the NS systems; instead, we briefly present the results of my previous research [9,10,11] on the three-dimensional NS system with various initial conditions in certain function spaces. Subsequently, we describe the methodology for investigating the studied n-dimensional NS problem with viscosity.

For this purpose, let us denote the velocity vector by $\mathsf{\nu }(x,t)\in R^n,(x\in R^n,t\ge 0)$ and the pressure by $P(x,t)\in R$. Then, the formulation of the n-dimensional incompressible NS system in vector form with the Cauchy initial condition in Cartesian coordinates is given by:

$${\displaystyle \frac{\partial \mathsf{\nu }}{\partial t}}+(\mathsf{\nu }\nabla )\mathsf{\nu }=f-{\displaystyle \frac{1}{\rho }}\nabla P+\mu \Delta \mathsf{\nu },$$

(1.1)

$$\mathrm{div}\mathsf{\nu }=0,$$

(1.2)

with initial conditions

$$\mathsf{\nu }{|}_{t=0}=\psi (x),\forall x\in R^n,$$

(1.3)

where $R^n∍\psi (x)$ is the known initial velocity vector, $R^n∍f(x,t)$ is external applied force (e.g. gravity), $(0,1)∍\mu$ is kinematic viscosity, $\rho$ is density, $\Delta$ is Laplace operator, $\nabla$ is Hamilton operator. These equations are to be solved for an unknown velocity vector $\mathsf{\nu }\in R^n$ and pressure $P(x,t)$, and equation (1.2) just says that the fluid is incompressible.

The aim of the present work is to establish the unique solvability of the original n-dimensional NS problem in a function space $G_{n,h}^1(D_0)$ equipped with the norm:

$$\left\{\begin{array}{l}{‖\mathsf{\nu }‖}_{G_{n,h}^1(D_0)}={\displaystyle \sum _{i=1}^n{‖{\mathsf{\nu }}_i‖}_{G_h^1(D_0)}=}{\displaystyle \sum _{i=1}^n\{}{\displaystyle \sum _{0\le \left|k\right|\le 2}^{\phantom{2}}{‖D_{\phantom{2}}^k{\mathsf{\nu }}_i‖}_{C(D_{+})}}+{‖{\mathsf{\nu }}_{it}‖}_{1h)}\},\\ {‖{\nu }_t‖}_{1h}=\underset{x\in R^n}{\sup }{\displaystyle \underset{0}{\overset{\infty }{\int }}h(s)}\left|{\nu }_s(x,s)\right|ds,(h(t)\in L^1(0,\infty );0\le h(t)\le h_1<\infty ,\forall t\in R_{+}),\\ {\displaystyle \underset{0}{\overset{\infty }{\int }}h(t)t^{m}dt\le h_2}<\infty ,(h_j=const,(j=1,2);h_0=\max (h_1,h_2);m=0,(-2^{-1})),\\ D_{+}=R^n\times R_{+},(R_{+}=[0,\infty ));D_0=R^n\times (0,\infty );D^0{\nu }_i={\nu }_i,(i=\overline{1,n}).\end{array}\right.$$

(1.4)

From the norm of the vector space defined in equation (1.4), it follows that the velocity vector of the NS problem (1.1)–(1.3) possesses second-order smoothness with respect to the spatial variables, and has a first-order partial derivative with respect to $t$, the variable $t>0$ provided that (a conditionally smooth solution) subject to the condition:

$${\displaystyle \underset{0}{\overset{\infty }{\int }}h(s)}\left|{\nu }_s(x,s)\right|ds<\infty ,\forall x\in R^n.$$

(1.5)

This implies that when the Cauchy condition is nonhomogeneous, as in (1.3), and taking into account that the NS system (1.1) is a nonlinear partial differential equation of parabolic type (or of a heat-conduction character; see, for example, Sobolev [12], Friedman and others), the introduction of a weighted vector space of the specified form $G_{n,h}^1(D_0)$ is a natural choice. The introduced space is not a Banach space (Sobolev [12], etc.).

Remark 1.  If condition (1.3) is replaced by $\mathsf{\nu }{|}_{t=0}=0,\forall x\in R^n,$ (1.3*),

then the velocity vector of the NSE (1.1) would possess second-order spatial smoothness and first-order smoothness with respect to the variable, i.e., it would be a smooth solution.

It should be noted that in my works [9,10], the 3D NSE for incompressible fluids were investigated under various values of the Cauchy condition and different time intervals. In these studies, different vector spaces were introduced depending on the specific value of the Cauchy condition. In particular, in [10], the 3D NSE was studied under condition (1.3*) in a Banach space $W_3(D_{*})$  equipped with the norm:

$$\left\{\begin{array}{l}{‖\mathsf{\nu }‖}_{W_3(D_{*})}={\displaystyle \sum _{i=1}^3{‖{\mathsf{\nu }}_i‖}_{W(D_{*})}=}{\displaystyle \sum _{i=1}^3\{}{\displaystyle \sum _{0\le \left|k\right|\le 2}^{\phantom{2}}{‖D_{\phantom{2}}^k{\mathsf{\nu }}_i‖}_{C(D_{*})}}+{‖{\mathsf{\nu }}_{it}‖}_{С(D_{*})}\},\\ D_{*}=R^3\times R_{+},(R_{+}=[0,\infty )).\end{array}\right.$$

(1.6)

In the article a 3D NS system with conditions of the form (1.2), (1.3) was investigated, where a method was proposed for constructing a solution (conditionally smooth in time) in a vector space $G_{3,h}^1(D_0=R^3\times (0,\infty ))$ with the norm

$$\left\{\begin{array}{l}{‖\mathsf{\nu }‖}_{G_{3,h}^1(D_0)}={\displaystyle \sum _{i=1}^3{‖{\mathsf{\nu }}_i‖}_{G_h^1(D_0)}=}{\displaystyle \sum _{i=1}^3\{}{\displaystyle \sum _{0\le \left|k\right|\le 2}^{\phantom{2}}{‖D_{\phantom{2}}^k{\mathsf{\nu }}_i‖}_{C(D)}}+{‖{\mathsf{\nu }}_{it}‖}_{1h)}\},(v\in R^3;\psi \in R^3),\\ {‖{\nu }_t‖}_{1h}=\underset{x\in R^3}{\sup }{\displaystyle \underset{0}{\overset{\infty }{\int }}h(s)}\left|{\nu }_s(x,s)\right|ds,(D=R^3\times R_{+};D_0=R^3\times (0,\infty );h(t)\in L^1(0,\infty )),\\ 0\le h(t)\le h_1=const<\infty ,\forall t\in R_{+}=[0,\infty ),\\ {\displaystyle \underset{0}{\overset{\infty }{\int }}h(s)s^{m}ds\le h_2}=const<\infty ,(h_0=\max (h_1,h_2);m=0,(-2^{-1})),\end{array}\right.$$

(1.7)

where $k=(k_1,k_2,k_3)$ is the multi-index,

$$\left\{\begin{array}{l}\mathsf{\nu }=({\mathsf{\nu }}_1,{\mathsf{\nu }}_2,{\mathsf{\nu }}_3),k=0:D^0{\mathsf{\nu }}_i\equiv {\mathsf{\nu }}_i,\\ k\ne 0:D^k{\mathsf{\nu }}_i={\displaystyle \frac{{\partial }^{\left|k\right|}{\mathsf{\nu }}_i}{\partial x_1{\phantom{2}}^{k_1}\partial x_2{\phantom{2}}^{k_2}\partial x_3{\phantom{2}}^{k_3}}},(i=\overline{1,3}),\left|k\right|={\displaystyle \sum _{j=1}^2{k}_j,}(k_j=0,1,2;j=1,2).\end{array}\right.$$

In this work, we follow the methodology of [11], but we immediately adapt the core idea of the method to the studied system of n-dimensional NSР (1.1) – (1.3), in order to avoid unnecessary expansion of the study.

2. Transformation of the NS System into Integral Form

In the theory of mathematical physics problems, various classes of partial differential equations have been studied in both bounded and unbounded domains under specific additional conditions (see Landau and Lifshith [4], Sobolev [12], Friedman [13], etc.]. For solving particular problems involving such equations, various mathematical transformations are employed that simplify the equations and enable one to find solutions in specific functional spaces — for instance, the Laplace and Fourier transforms, Riemann and Green functions, and others. However, in the general case of the NS system with conditions (1.2) and (1.3), these tools are not applicable.

Therefore, in this section, the proposed method for analyzing the n-dimensional NS system for incompressible fluids with low viscosity — without invoking any additional conditions — transforms the studied system into an integral form.

To this end, let the velocity vector at the initial moment of time, $\psi (x)\in R^n$, and the external force, $f(x,t)\in R^n$, satisfy the following conditions:

$$\left\{\begin{array}{l}{\displaystyle \underset{R^n}{\int }\psi }(\tau )d\tau +{\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{R^n}{\int }f}}(\tau ,s)d\tau ds=\lambda ,(\tau \in R^n),\\ \psi (x)\in C^2(R^n),f(x,t)\in C^{1,0}(D_{+}),(x\in R^n;{С}^2(R^n)\equiv C^{2,\mathrm{...},2}(R^n)),\end{array}\right.$$

(2.1)

here $R^n∍\lambda$ is a known vector with positive constant components: $0<{\lambda }_i,(i=\overline{1,n})$. Then we can find the velocity vector $v(x,t)\in R^n$ based on the transformation:

$$\mathsf{\nu }=\theta \lambda +(\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}}))J(x,t),$$

(2.2)

where $R^n∍J$ is the given vector:

$$\left\{\begin{array}{l}J(x,t)\equiv {\displaystyle \frac{1}{2^n{(\sqrt{\mu \pi t})}^n}}{\displaystyle \underset{R^n}{\int }\psi (\tau )\exp }(-{\displaystyle \frac{{\left|x-\tau \right|}^2}{4\mu t}})d\tau ={\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{R^n}{\int }\psi (x+2\xi \sqrt{\mu t})\times }\\ \times \exp (-{\left|\xi \right|}^2)d\xi ,(x,\tau ,\xi \in R^n),\\ \left|x-\tau \right|=\sqrt{{\displaystyle \sum _{i=1}^n{(x_i-{\tau }_i)}^2}};(0<\mu <1;0<{\delta }_0=const<1),\\ {\left.J\right|}_{t=0}=\psi (x),\forall x\in R^n,\\ G(x,\tau ,t)\equiv {\displaystyle \frac{1}{2^n{(\sqrt{\mu \pi t})}^n}}\exp (-{\displaystyle \frac{{\left|x-\tau \right|}^2}{4\mu t}}),(t>0),\\ L[G]\equiv {\displaystyle \frac{\partial G}{\partial t}}-\mu \Delta G=0,\end{array}\right.$$

(2.3)

in this case, $\theta (x,t)-$ a new unknown scalar function with the condition:

$$\theta (x,0)=0,\forall x\in R^n,$$

(2.4)

$0<{\delta }_0$ is the introduced known constant, ensuring the application of the Picard method for the obtained system IE-2, which follows from the original problem taking into account (2.2).

In addition, taking into account condition (1.2) from (2.2) we obtain

$$\left\{\begin{array}{l}\mathrm{div}\mathsf{\nu }=0,(\mathrm{div}\psi =0):\mathrm{div}J={\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{R^n}{\int }\exp }(-{\left|\xi \right|}^2)\mathrm{div}\psi (x+2\xi \sqrt{\mu t})d\xi =0,\\ 0=\mathrm{div}\mathsf{\nu }=\mathrm{div}(\theta \lambda )+(\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}}))\mathrm{div}J=\mathrm{div}(\theta \lambda )={\displaystyle \sum _{i=1}^n{\theta }_{x_i}}{\lambda }_i.\end{array}\right.$$

(2.5)

Hence, from conditions (2.3), (2.4), and (2.5), it follows that the introduced algorithm (2.2) is consistent with conditions (1.2) and (1.3).

We will now prove the following results:

a.1) Under conditions (2.3), (2.4), and (2.5), based on the transformation (2.2), the nonlinear convective terms in the NS system are linearized;

a.2) The algorithm (2.2) transforms the nonlinear NS system into a linear form, from which the Poisson equation for pressure is derived;

a.3) Eliminating the pressure from the resulting linear system leads to an inhomogeneous linear heat-type equation with condition (2.4) for the function $\theta (x,t)$;

a.4) The Cauchy problem obtained in point (a.3) is transformed into a system of integral equations of the second kind (IE-2), where the integrals over the spatial coordinates belong to the class of Poisson integrals. From the unique solvability of this system, the unique solvability of the original problem in $G_{n,h}^1(D_0)$ follows.

Note 1. All the results in points (a.1–a.4) will be rigorously justified later. It should be noted that the first condition (a.1) is novel in the theory of the NS system for incompressible fluids with low viscosity (except for my own works [9, 10, and 11] on the 3D NSE). In our case, the linearization of the convective terms is not associated with the assumptions made in works such as [4 and 6], etc. In those studies, the linearization of the NSE was achieved by neglecting the nonlinear convective terms based on the assumption that these terms are sufficiently small compared to viscous forces, which is valid when the viscosity is relatively large, and so on.

In my research, the linearization of the convective terms is performed based on (2.2), taking into account condition (1.2), and this condition is a key factor in the theory of NS systems. Many researchers who study my work conclude that this result is either impossible or artificial — but that is not the case. This result is not applicable to arbitrary nonlinear differential equations. However, the equation under study in this paper is the NSE with condition (1.2), and this condition — in combination with (2.2) — allows for the linearization of the convective terms, and hence, the entire NS system. (This Note 1 is provided to clarify the essence of the present work.)

2.1. Linearization of the Convective Terms in the NSE

To prove statement (a.1), we first formulate the following lemma:

Lemma 1. In the case of (2.2), when the incompressibility condition (1.2) holds, the nonlinear convective terms in the NS system (1.1) are linearized with respect to the functions $\theta ,{\theta }_{x_i},(i=\overline{1,n})$.

Proof. Since, under condition (1.2), the algorithm (2.2) satisfies conditions (2.3), (2.4), and (2.5), substituting (2.2) into the NS system (1.1) for the inertial terms yields:

$$(\mathsf{\nu }\nabla )\mathsf{\nu }=(\theta \lambda \nabla )\theta \lambda +(\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}}))[(\theta \lambda \nabla )J+(J\nabla )\theta \lambda ]+(\exp (-{\displaystyle \frac{2t}{\mu {\delta }_0}}))(J\nabla )J.$$

(2.6)

From this it is clear that if condition (2.5) is met, the equality follows

$$(\theta \lambda \nabla )\theta \lambda ={\lambda }_i\theta ({\displaystyle \sum _{j=1}^n{\lambda }_j}{\theta }_{x_j})=0,(i=\overline{1,n}).$$

(2.7)

Therefore, from (2.6) we obtain

$$(\mathsf{\nu }\nabla )\mathsf{\nu }=(\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}}))[(\theta \lambda \nabla )J+(J\nabla )\theta \lambda ]+(\exp (-{\displaystyle \frac{2t}{\mu {\delta }_0}}))(J\nabla )J.$$

(2.8)

This means that the inertial terms of the NS system (1.1), indeed — taking into account (1.2) and (2.2) — are linearized with respect to the newly introduced function and its partial derivatives with respect to $x\in R^n$, while the nonlinearity is transferred to the known vector function $J(x,t)$ and its partial derivatives with respect to $x\in R^n.$ Thus, Lemma 1 is proved.

2.2. Proof of Condition a.2)

To prove condition (a.2), taking into account the assumptions of Lemma 1 and substituting (2.2) into the NS system (1.1), we obtain a system of linear partial differential equations of the heat conduction type with variable coefficients:

$$\begin{array}{l}{\displaystyle \frac{\partial \theta }{\partial t}}\lambda +(\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}}))[(\theta \lambda \nabla )J+(J\nabla )\theta \lambda ]={\displaystyle \frac{1}{\mu {\delta }_0}}(\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}}))J-(\exp (-{\displaystyle \frac{2t}{\mu {\delta }_0}}))\times \\ \times (J\nabla )J+f-{\rho }^{-1}\nabla P+(\mu \Delta \theta )\lambda .\end{array}$$

(2.9)

Then, applying the operation $\mathrm{div}$ to system (2.9), we derive the Poisson equation for pressure:

$$\left\{\begin{array}{l}{\displaystyle \frac{1}{\rho }}\Delta P=-{\displaystyle \sum _{i=1}^n{\displaystyle \sum _{k=1}^n{\mathsf{\nu }}_{i{x}_k}}}{\mathsf{\nu }}_{k{x}_i}=-{\chi }_n^{\phantom{2}}\{F_0+{\chi }_n^{-1}(\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}}))[{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{\lambda }_j{J}_{i{x}_j}){\theta }_{x_i}}}+\\ +{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{J}_{j{x}_i}{\theta }_{x_j}}}){\lambda }_i]\}\equiv -{\chi }_n\Omega (x,t),\\ {\chi }_n=2(n-2)(\Gamma (2^{-1}n){)}^{-1}\sqrt{{\pi }^n},(n\ge 3),\end{array}\right.$$

(2.10)

since the following take place:

$$\left\{\begin{array}{l}\mathrm{div}f=0;\mathrm{div}J=0;\mathrm{div}({\theta }_t\lambda )=0;\mathrm{div}(\mu \Delta \theta )\lambda =0;F_0\equiv {\chi }_n^{-1}\exp (-{\displaystyle \frac{2t}{{\delta }_0\mu }}){\displaystyle \sum _{i=1}^n{\displaystyle \sum _{j=1}^n{J}_{i{x}_j}{J}_{j{x}_i},}}\\ \mathrm{div}\{\exp (-{\displaystyle \frac{t}{{\delta }_0\mu }})[(\theta \lambda \nabla )J+(J\nabla )\theta \lambda ]+\exp (-{\displaystyle \frac{2t}{{\delta }_0\mu }})(J\nabla )J\}=F_0+\exp (-{\displaystyle \frac{t}{{\delta }_0\mu }})\times \\ \times ({\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{\lambda }_j{J}_{i{x}_j}){\theta }_{x_i}}}+{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{J}_{j{x}_i}{\theta }_{x_j}}}){\lambda }_i).\end{array}\right.$$

Further, from DE (2.10) it follows:

$$\left\{\begin{array}{l}\Omega (x,t)\equiv F_0+{\chi }_n^{-1}\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}})[{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{\lambda }_j{J}_{i{x}_j}){\theta }_{x_i}}}+{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{J}_{j{x}_i}{\theta }_{x_j}}}){\lambda }_i],\\ P(x,t)=\rho {\displaystyle \underset{R^n}{\int }\frac{1}{r^{n-2}}\Omega }(\tau ,t)d\tau ,(x,\tau \in R^n,r=\left|x-\tau \right|),\end{array}\right.$$

(2.11)

in this case the following takes place:

$${\displaystyle \frac{\partial }{\partial x}}P=\rho {\displaystyle \underset{R^n}{\int }\Omega (\tau ,t)\frac{(\tau -x)(n-2)}{r^n}d\tau },(\tau -x\in R^n).$$

(2.12)

Equation (2.11) is referred to as the Newtonian potential (see Sobolev [12], еtс.). Moreover, the solution to the Poisson equation (2.10) that tends to zero at infinity will be unique, provided that the functions ${\theta }_{x_i},(i=\overline{1,n})$ are unique, since contains $\Omega$ these functions. Thus, condition (a.2) is proved.

2.3. Proof of Condition a.3)

To prove condition (a.3), we first examine system (2.9), since it contains the scalar function $\theta (x,t)$.

Remark 2. To better understand (2.9), we present a specific example from the field of systems of linear algebraic equations (SLAE), in which the system contains a single unknown variable z (for simplicity, we consider a square system ($3\times 3$)), i.e.:

${\lambda }_{i}z+a_i={\lambda }_{i}b,(0<{\lambda }_i,i=\overline{1,3}).$ (0.1)

It follows from this that, if the following holds

${\lambda }_1^{-1}{a}_1={\lambda }_2^{-1}{a}_2={\lambda }_3^{-1}{a}_3=a_0$, (0.2)

then z is uniquely determined in the form

$z=b-a_0.$(0.3)

But we can define z differently, i.e.:

$z=b-{({\lambda }_1+{\lambda }_2+{\lambda }_3)}^{-1}(a_1+a_2+a_3)$ (0.4)

Alternatively, with respect to (0.4), performing a certain mathematical transformation, we obtain

$z=b-{({\lambda }_1+{\lambda }_2+{\lambda }_3)}^{-1}[{\lambda }_1({\lambda }_1^{-1}{a}_1)+{\lambda }_2({\lambda }_2^{-1}{a}_2)+{\lambda }_3({\lambda }_3^{-1}{a}_3)]\underset{(0.2)}{=}b-a_0.$ (0.5)

This means that the first and second paths are equivalent, i.e. under condition (0.2), z is in fact uniquely determined from (0.1)

Therefore, using the idea of ​​system (0.1) relative to system (2.9), the condition is allowed:

$$\left\{\begin{array}{l}{\lambda }_1^{-1}{J}_1\equiv {\lambda }_2^{-1}{J}_2\equiv \mathrm{...}\equiv {\lambda }_n^{-1}{J}_n;{\lambda }_1^{-1}({\displaystyle \sum _{j=1}^n{\lambda }_j}{J}_{1x_j})\equiv {\lambda }_2^{-1}({\displaystyle \sum _{j=1}^n{\lambda }_j{J}_{2x_j}})\equiv \mathrm{...}\equiv ({\lambda }_n^{-1}{\displaystyle \sum _{j=1}^n\phantom{2}}{\lambda }_j{J}_{n{x}_j}),\\ {\lambda }_1^{-1}\{{\displaystyle \frac{1}{\rho }}{P}_{x_1}-f_1+(\exp (-{\displaystyle \frac{2t}{\mu {\delta }_0}})){\displaystyle \sum _{j=1}^n{J}_j{J}_{1x_j}}\}\equiv \mathrm{...}\equiv {\lambda }_n^{-1}\{{\displaystyle \frac{1}{\rho }}{P}_{x_n}-f_n+(\exp (-{\displaystyle \frac{2t}{\mu {\delta }_0}})){\displaystyle \sum _{j=1}^n{J}_j}{J}_{n{x}_j}\}.\end{array}\right.$$

(2.13)

The final condition (2.13) represents a condition of unique compatibility for system (2.9), since $\theta (x,t)$ is a scalar function. It should be noted that (2.13) is a condition of the type (0.2) for the studied system with the scalar function $\theta (x,t)$. Moreover, when solving system (2.9), we follow the second approach, as demonstrated in the case of (0.4).

Indeed, substituting (2.12) into system (2.9) and taking into account: $d_0={\lambda }_1+\mathrm{...}+{\lambda }_n\ne 0$, we obtain a Cauchy problem for an inhomogeneous linear heat-type equation with variable coefficients of the form:

$$\left\{\begin{array}{l}{\theta }_t=\Phi +\zeta (x,t)\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}})+\mu \Delta \theta ,\\ \zeta (x,t)=-\{d_0^{-1}\theta (.){\displaystyle \sum _{i=1}^n(}{\displaystyle \sum _{j=1}^n{\lambda }_j}{I}_{i{x}_j}(.))+{\displaystyle \sum _{j=1}^n{\theta }_{x_j}}(.)I_j(.)+d_0^{-1}{\chi }_n^{-1}({\displaystyle \underset{R^n}{\int }{\displaystyle \sum _{k=1}^n\frac{(n-2)}{r_1^n}{\overline{\xi }}_k}}\times \\ \times [{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{\lambda }_j{I}_{i{\tau }_j}(x}}+\overline{\xi },t)){\theta }_{{\tau }_i}(x+\overline{\xi },t)+{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{I}_{j{\tau }_i}(x}}+\overline{\xi },t){\theta }_{{\tau }_j}(x+\overline{\xi },t)){\lambda }_i]d\overline{\xi })\},\\ \theta (x,t)\left|{\phantom{2}}_{t=0}\right.=0,\forall x\in R^n,\end{array}\right.$$

(2.14)

where

$$\left\{\begin{array}{l}\Phi \equiv {\displaystyle \sum _{i=1}^n{\Phi }_i};{\Phi }_1\equiv d_0^{-1}\left\{\phantom{2}\right.{\displaystyle \sum _{i=1}^n{f}_i}-(\exp (-{\displaystyle \frac{2t}{\mu {\delta }_0}})){\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n\phantom{2}}}{J}_j{J}_{i{x}_j})\left.\phantom{2}\right\},\\ {\Phi }_2\equiv d_0^{-1}\left\{{\displaystyle \frac{1}{\mu {\delta }_0}}\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}})\left.{\displaystyle \sum _{i=1}^n{J}_i}\right\}\right.,\\ {\Phi }_3\equiv d_0^{-1}[-{\displaystyle \underset{R^n}{\int }{\displaystyle \sum _{i=1}^n\frac{{\overline{\xi }}_i(n-2)}{r_1^n}{F}_0}}(x+\overline{\xi };t)d\overline{\xi }]=d_0^{-1}(\exp (-{\displaystyle \frac{2t}{\mu {\delta }_0}}\left)\right)[-{\chi }_n^{-1}\times \\ \times {\displaystyle \underset{R^n}{\int }{\displaystyle \sum _{i=1}^n\frac{{\overline{\xi }}_i(n-2)}{r_1^n}}}{\displaystyle \sum _{m=1}^n{\displaystyle \sum _{k=1}^n(J_{m{\tau }_k}}}(x+\overline{\xi };t)J_{k{\tau }_m}(x+\overline{\xi };t))d\overline{\xi }],(r_1=({\displaystyle \sum _{n=1}^n{\overline{\xi }}_n^2{)}^{2^{-1}}).}\end{array}\right.$$

(2.15)

Condition (a.3) is proved.

2.4. Proof of Condition a.4)

It is known that the Cauchy problem (2.14) with sufficiently smooth initial data is uniquely solvable [12 and 13]; in our case, in the space $G_h^1(D_0)$ . To prove this, we first transform the Cauchy problem (2.14) into an IE-2 system, from which the unique solvability of this IE system implies the unique solvability of the original NS problem in $G_{n,h}^1(D_0)$, taking into account (2.2). Indeed, since the Cauchy problem (2.14) implies:

$$\left\{\begin{array}{l}\theta ={\rm Y}+{\displaystyle \frac{1}{2^n\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\int }(\exp }}(-{\displaystyle \frac{r^2}{4\mu (t-s)}}))(\exp (-{\displaystyle \frac{s}{\mu {\delta }_0}}))\zeta (\tau ,s){\displaystyle \frac{d\tau ds}{{(\sqrt{\mu (t-s)})}^n}}=\\ ={\rm Y}+{\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }}(-({\left|\xi \right|}^2+{\displaystyle \frac{s}{\mu {\delta }_0}})))\zeta (x+2\xi \sqrt{\mu (t-s)},s)d\xi ds\equiv ({\Gamma }_0\zeta )(x,t),\\ {\theta }_{x_i}={{\rm Y}}_{x_i}+{\displaystyle \frac{1}{2^n\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\int }(\exp }}(-({\displaystyle \frac{r^2}{4\mu (t-s)}}+{\displaystyle \frac{s}{\mu {\delta }_0}}))){\displaystyle \frac{-(x_i-{\tau }_i)}{2\mu (t-s)}}\zeta (\tau ,s)\times \\ \times {\displaystyle \frac{d\tau ds}{{(\sqrt{\mu (t-s)})}^n}}={{\rm Y}}_{x_i}+{\displaystyle \frac{1}{\sqrt{{\pi }^3}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }}(-({\left|\xi \right|}^2+{\displaystyle \frac{s}{\mu {\delta }_0}})))\zeta (x+2\xi \sqrt{\mu (t-s)},s)\times \\ \times {\displaystyle \frac{{\xi }_{i}d\xi ds}{\sqrt{\mu (t-s})}}\equiv ({\Gamma }_i\zeta )(x,t),(i=\overline{1,n}),\\ \zeta (x,t)\in C^{1,0}(D_{+}),(C^{1,0}(D_{+})\equiv C^{1,\mathrm{...},1,0}(D_{+}),x\in R^n);{\rm Y}\equiv {\displaystyle \sum _{i=1}^2\phantom{2}}{{\rm Y}}_i,\\ {{\rm Y}}_1={\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }}(-{\left|\xi \right|}^2))({\Phi }_1(x+2\xi \sqrt{\mu (t-s)},s)+{\Phi }_3(x+2\xi \sqrt{\mu (t-s)},s))d\xi ds,\\ {{\rm Y}}_2={\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\int }(\exp }}(-{\left|\xi \right|}^2)){\Phi }_2(x+2\xi \sqrt{\mu (t-s)},s)d\xi ds,\end{array}\right.$$

(2.16)

then, taking into account the auxiliary function $\zeta (x,t)$, we obtain the system:

$$\left\{\begin{array}{l}\theta =({\Gamma }_0\zeta )(x,t),\\ \zeta =-\{d_0^{-1}({\Gamma }_0\zeta ){\displaystyle \sum _{i=1}^n(}{\displaystyle \sum _{j=1}^n{\lambda }_j}{I}_{i{x}_j}(.))+{\displaystyle \sum _{j=1}^n({\Gamma }_j\zeta )}{I}_j(.)+d_0^{-1}({\chi }_n^{-1}{\displaystyle \underset{R^n}{\int }{\displaystyle \sum _{k=1}^n\frac{{\tilde{\xi }}_k(n-2)}{r_1^n}\times }}\\ \times [{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{\lambda }_j{I}_{i{\tau }_j}(x}}+\tilde{\xi },t))({\Gamma }_i\zeta )(x+\tilde{\xi },t)+{\displaystyle \sum _{i=1}^n({\displaystyle \sum _{j=1}^n{I}_{j{\tau }_i}(x}}+\tilde{\xi },t)({\Gamma }_j\zeta )(x+\\ +\tilde{\xi },t)){\lambda }_i]d\tilde{\xi })\}\equiv (\Gamma \zeta )(x,t),(\tau =x+\tilde{\xi }\in R^n),\end{array}\right.$$

(2.17)

in this case, the following conditions are met with respect to the known functions of the system (2.17):

$$\left\{\begin{array}{l}\left|D^k{{\rm Y}}_2\right|\le {\beta }_1,\left|D^{k}J\right|\le {\beta }_2,\forall (x,t)\in D_{+},\\ {{\rm Y}}_{2x_i}={\displaystyle \frac{d_0^{-1}}{\mu {\delta }_0\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }}(-({\left|\xi \right|}^2+{\displaystyle \frac{s}{\mu {\delta }_0}})))({\displaystyle \sum _{j=1}^n{J}_{j{l}_i}}(x+2\xi \sqrt{\mu (t-s)},s))d\xi ds\equiv {\Phi }_{4,i},\\ {{\rm Y}}_{2x_i^2}={\displaystyle \frac{\partial }{\partial x_i}}({\Phi }_{4,i}(x,t)),(i=\overline{1,n};l=x+2\xi \sqrt{\mu (t-s)}\in R^n),\end{array}\right.$$

$$\left\{\begin{array}{l}{\displaystyle \frac{1}{\mu {\delta }_0}}{\displaystyle \underset{0}{\overset{t}{\int }}\exp (-\frac{s}{\mu {\delta }_0})ds=1-\exp (-\frac{t}{\mu {\delta }_0})\le 1,\forall t\in R_{+}},\\ {{\rm Y}}_{2t}={\Phi }_2+{\displaystyle \frac{1}{\mu {\delta }_0\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }}(-({\left|\xi \right|}^2+{\displaystyle \frac{s}{\mu {\delta }_0}}))){\displaystyle \sum _{k=1}^n\frac{{\xi }_k\sqrt{\mu }}{\sqrt{t-s}}}{J}_{l_k}(x+2\xi \sqrt{\mu (t-s)},s)d\xi ds,\\ {‖{\Phi }_2‖}_{1h}=\underset{R^n}{\sup }{\displaystyle \underset{0}{\overset{\infty }{\int }}h(s)\left|{\Phi }_2(x,s)\right|}ds\le d_0^{-1}{\beta }_2{h}_0={\beta }_3,\\ {\displaystyle \frac{1}{\sqrt[\phantom{2}]{\mu }}}{\displaystyle \underset{0}{\overset{t}{\int }}(\exp (-\frac{s}{\mu {\delta }_0}))\frac{ds}{\sqrt{t-s}}\le }{\displaystyle \frac{1}{\sqrt{\mu }}}({\displaystyle \underset{0}{\overset{t}{\int }}(\exp (-\frac{2(\sqrt{t}-\sqrt{\tau })(\sqrt{t}+\sqrt{\tau })}{\mu {\delta }_0}}){\displaystyle \frac{d\tau }{\sqrt{\tau }}}{)}^{{\displaystyle \frac{1}{2}}}({\displaystyle \underset{0}{\overset{t}{\int }}\frac{d\tau }{\sqrt{\tau }}{)}^{\frac{1}{2}}}\le \\ \le \sqrt{2{\delta }_0}[{\displaystyle \underset{0}{\overset{t}{\int }}(\exp (-\frac{2(\sqrt{t}-\sqrt{\tau })\sqrt{t}}{\mu {\delta }_0}))\frac{\sqrt{t}d\tau }{\mu {\delta }_0\sqrt{\tau }}{]}^{\frac{1}{2}}=}\sqrt{2{\delta }_0}[{\displaystyle \underset{0}{\overset{t}{\int }}(\exp }(-{\displaystyle \frac{2(\sqrt{t}-\sqrt{\tau })\sqrt{t})}{\mu {\delta }_0}})\times \\ \times d(-{\displaystyle \frac{2}{\mu {\delta }_0}}(\sqrt{t}-\sqrt{\tau })\sqrt{t}){]}^{{\displaystyle \frac{1}{2}}}\le \sqrt{2{\delta }_0}{[{\displaystyle \underset{0}{\overset{\infty }{\int }}\phantom{2}}\exp (-\rho )d\rho ]}^{{\displaystyle \frac{1}{2}}}\le \sqrt{2{\delta }_0},\\ \underset{{\mathrm{D}}_{+}}{\sup }{\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }(-{\left|\xi \right|}^2)){\displaystyle \sum _{k=1}^n\left|{\xi }_k\right|\times \left|J_{l_k}(x+2\xi \sqrt{\mu (t-s)},s)\right|d\xi }\le {\beta }_2{\beta }_4={\beta }_5,\\ \underset{{\mathrm{D}}_{+}}{\sup }{\displaystyle \frac{1}{\sqrt{\mu {\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}\frac{1}{\sqrt{t-s}}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }}(-({\left|\xi \right|}^2+{\displaystyle \frac{s}{\mu {\delta }_0}}))){\displaystyle \sum _{k=1}^n\left|{\xi }_k\right|\times \left|J_{l_k}\right|d\xi ds}\le {\beta }_5\sqrt{2{\delta }_0}={\beta }_6,\\ {‖{{\rm Y}}_{2t}‖}_{1h}\le {\beta }_3+{\beta }_6\le {\beta }_7,\\ {‖{{\rm Y}}_2‖}_{G_h^1(D_0)}={\displaystyle \sum _{0\le \left|k\le 2\right|}{‖{{\rm Y}}_2‖}_{C(D_{+})}+{‖{{\rm Y}}_{2t}‖}_{1h}}\le {\beta }_8,(0<{\beta }_i=const;i=\overline{1,8}).\end{array}\right.$$

(2.18)

Similar assessments can be made regarding ${{\rm Y}}_1$ in $G_h^1(D_0)$, therefore we have ${\rm Y}\in G_h^1(D_0).$

Furthermore, estimates (2.18) provide methods for assessing the regularity of certain integrals with respect to small viscosity, which may appear in the estimates of system (2.17). Therefore, taking (2.18) into account, we can assert that the linear integral operators ${\Gamma }_0,{\Gamma }_i,(i=\overline{1,n})$ introduced in system (2.17) are regular with respect to small viscosity $\mu$. Then, for the second integral equation of the second kind (IE-2) in system (2.17), the conditions of Banach’s fixed-point principle are satisfied, since in the corresponding estimates, the Lipschitz constant of the introduced linear integral operator $\Gamma$ involves known small quantities ${\delta }_0$ and $\sqrt{{\delta }_0}$, i.e.:

$$\left\{\begin{array}{l}\left|\Gamma \zeta -\Gamma {\zeta }_0\right|\le L_{\Gamma }{‖\zeta -{\zeta }_0‖}_{C(D_{+})},\\ L_{Г}=k_1{\delta }_0+k_2\sqrt{{\delta }_0}\le \overline{k}\sqrt{{\delta }_0}<1,(2^{-1}\overline{k}=\max (k_1,k_2)),\\ 1)\overline{k}\le 1,0<{\delta }_0=const<1\mathrm{or}2)1<\overline{k}\le k_0,0<{\delta }_0<k_0{\phantom{2}}^{-2}<1.\end{array}\right.$$

(2.19)

This means that the given operator is a contraction operator; therefore, we can formulate the following lemma:

Lemma 2. Under conditions (2.18), (2.19) and

$$\left\{\begin{array}{l}{‖\Gamma {\zeta }_0-{\zeta }_0‖}_{C(D_{+})}\le (1-L_{\Gamma })r_0,\\ S_{r_0}({\zeta }_0)=\{\zeta :\left|\zeta -{\zeta }_0\right|\le r_0,\forall (x,t)\in D_{+}\}.\end{array}\right.$$

(2.20)

The second integral equation (IE) of system (2.17) is uniquely solvable in $C^{1,0}(D_{+})$, and its solution is constructed according to Picard's iteration method:

$${\zeta }_{m+1}=\Gamma {\zeta }_m,(m=0,1,\mathrm{...})$$

(2.21)

with an estimate

$${‖{\zeta }_{m+1}-\zeta ‖}_{C(D_{+})}\le {(L_{Г})}^{m+1}{r}_0\underset{L_{Г}<1,(m\to \infty )}{\to }0$$

(2.22)

Therefore, the function $\theta (x,t)$ is also uniquely determined as the solution of the first integral equation of system (2.17) in $G_h^1(D_0).$

Proof. (A) Since, by the assumption of Lemma 2, condition (2.20) holds, it follows that the following inequality is satisfied:

$$\begin{array}{l}{‖\Gamma \zeta -{\zeta }_0‖}_{C(D_{+})}\le {‖\Gamma \zeta -\Gamma {\zeta }_0‖}_{C(D_{+})}+{‖\Gamma {\zeta }_0-{\zeta }_0‖}_{C(D_{+})}\le L_{\Gamma }{‖\zeta -{\zeta }_0‖}_{C(D_{+})}+(1-L_{{Г}_0})r_0\le \\ \le L_{\Gamma }{r}_0+(1-L_{\Gamma })r_0=r_0,\end{array}$$

this implies that $\Gamma :S_{r_0}({\zeta }_0)\to S_{r_0}({\zeta }_0),$ (2.23)

i.e., the operator Г maps its domain into itself. Therefore, the conditions of the Banach fixed-point principle are satisfied for the operator (see (2.19), (2.23)). It follows that the integral equation involving admits a solution, and moreover, this solution is unique in $C^{1,0}(D_{+})$ (the first-order partial derivatives with respect to the spatial variables are continuous for $\forall (x,t)\in D_{+})$.

Uniqueness is proved by contradiction. Indeed, suppose that the system of equations (2.17) admits not only the solution $\zeta$, but also another solution $\overline{\zeta }$, where

$\overline{\zeta }=(\Gamma \overline{\zeta })(x,t).$

It follows:

$\left|\zeta -\overline{\zeta }\right|=\left|\Gamma \zeta -\Gamma \overline{\zeta }\right|\le L_{\Gamma }{‖\zeta -\overline{\zeta }‖}_C$or

${‖\zeta -\overline{\zeta }‖}_{C(D_{+})}\le L_{\Gamma }{‖\zeta -\overline{\zeta }‖}_{C(D_{+})},(L_{\Gamma }<1),$

which is satisfied if and only if:

${‖\zeta -\overline{\zeta }‖}_{C(D_{+})}=0,$

i.e. $\zeta \equiv \overline{\zeta },$ and this guarantees the uniqueness of the function $\zeta (x,t)$. Which was required to show.

Furthermore, the solution of the integral equation with operator Г, as noted above, is obtained according to Picard's iteration scheme (2.21). Taking into account the conclusions derived from this method [12], we obtain:

$$\left\{\begin{array}{l}{‖{\zeta }_{m+1}-{\zeta }_m‖}_{C(D_{+})}\le L_{\Gamma }{‖{\zeta }_m-{\zeta }_{m-1}‖}_{C(D_{+})}\le \mathrm{...}\le {(L_{\Gamma })}^m{‖{\zeta }_1-{\zeta }_0‖}_{C(D_{+})}\le {(L_{\Gamma })}^m{r}_{*},\\ {‖{\zeta }_{m+k}-{\zeta }_m‖}_{C(D_{+})}\le L_{\Gamma }{\displaystyle \sum _{j=0}^{k-1}\phantom{2}}{‖{\zeta }_{m+j}-{\zeta }_{m+j-1}‖}_{C(D_{+})}\le \mathrm{...}\le {(L_{\Gamma })}^m{(1-L_{\Gamma })}^{-1}{r}_{*},\\ \left|{\zeta }_1-{\zeta }_0\right|\le r_{*},\forall (x,t)\in D_{+},\end{array}\right.$$

(2.24)

where ${\zeta }_0$ is initial estimate.

Means, based on the conclusions of the Picard's method we obtain:

$$\left\{\begin{array}{l}{‖\zeta -{\zeta }_m‖}_{C(D_{+})}\le {(L_{\Gamma })}^m{r}_0\underset{L_{\Gamma }<1,(m\to \infty )}{\to }0,\\ {‖\zeta ‖}_{C(D_{+})}\le {(1-L_{\Gamma })}^{-1}{M}_1=M_2,\\ {‖\zeta ‖}_{C^{1,0}(D_{+})}\le M_3,\end{array}\right.$$

(2.25)

i.e. the second equation (SE) of the system (2.17) is solvable in $C^{1,0}(D_{+}).$ The first part of Lemma 2 has been established.

(B) We now proceed to prove the second part of Lemma 2. Using system (2.17), we show that the function $\theta (x,t)$ is uniquely determined by the estimate.

$$\left\{\begin{array}{l}{\theta }_m=({\Gamma }_0{\zeta }_m)(x,t),\\ {\zeta }_{m+1}=(Г{\zeta }_m)(x,t),(m=0,1,2,\mathrm{...}),\\ {‖{\theta }_m-\theta ‖}_{C(D_{+})}\le \mu {\delta }_0{‖{\zeta }_m-\zeta ‖}_{C(D_{+})}\le \mu {\delta }_0{L}_{{\phantom{2}}^{Г}}^m{r}_0\underset{L_{Г}<1,m\to \infty }{\to }0,\\ {‖{\zeta }_m-\zeta ‖}_{C(D_{+})}\le L_{{\phantom{2}}^{Г}}^m{r}_0,\\ L_{Г}=\overline{k}\sqrt{{\delta }_0}<1,(\left|\zeta \right|\le M_2,\forall (x,t)\in D_{+}),\\ {‖\theta ‖}_{C(D_{+})}\le {‖{\rm Y}‖}_{C(D_{+})}+\mu {\delta }_0{M}_2\le {\beta }_9.\end{array}\right.$$

(2.26)

Further, since equations (2.16), (2.17) contain:

$$\theta (x,t),{\theta }_{x_i}(x,t),(i=\overline{1,n}),\zeta (x,t)\in C^{1,0}(D_{+}),$$

then, taking into account

$$\left\{\begin{array}{l}{\theta }_{x_i^2}={{\rm Y}}_{x_i^2}+{\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }}(-({\left|\xi \right|}^2+{\displaystyle \frac{s}{\mu {\delta }_0}}))){\zeta }_{l_i}(x+2\xi \sqrt{\mu (t-s)},s){\displaystyle \frac{{\xi }_i}{\sqrt{\mu (t-s)}}}d\xi ds,\\ {\theta }_t={{\rm Y}}_t+\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}})\zeta (x,t)+{\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{0}{\overset{t}{\int }}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }}(-({\left|\xi \right|}^2+{\displaystyle \frac{s}{\mu {\delta }_0}})))\sqrt{\mu }{\displaystyle \sum _{j=1}^n\frac{{\xi }_j}{\sqrt{t-s}}}{\zeta }_{l_j}(x+\\ +2\xi \sqrt{\mu (t-s)},s)d\xi ds,(i=\overline{1,n};l=x+2\xi \sqrt{\mu (t-s)}\in R^n),\end{array}\right.$$

(2.27)

moreover, from the estimate systems (2.16), (2.27) on the basis of condition (2.18) it follows:

$${‖\theta ‖}_{G_h^1(D_0)}={\displaystyle \sum _{0\le \left|k\right|\le 2}^{\phantom{2}}{‖D_{\phantom{2}}^k\theta ‖}_{C(D_{+})}}+{‖{\theta }_t(x,t)‖}_{1h}\le {\beta }_{10}=const.$$

(2.28)

Lemma 2 is now proved.

We note that in step (a.2), applying $\mathrm{div}$ to system (2.9) yielded a Poisson equation for the pressure. This equation was then transformed into the form (2.11), where the function $\Omega (x,t)$ involves the functions ${\theta }_{x_i},(i=\overline{1,n})$. Therefore, based on Lemma 2, we conclude that in this case the pressure is determined, since the right-hand side of equation (2.11) is a known function. This completes the proof of the claim.

Finally, taking into account the transformation (2.2), we obtain

$$\left\{\begin{array}{l}{\mathsf{\nu }}_{i,m}={\theta }_m{\lambda }_i+\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}})J_i(x,t),(i=\overline{1,n};m=0,1,2,\mathrm{...}),\\ {‖{\mathsf{\nu }}_{i,m}-{\mathsf{\nu }}_i‖}_{C(D_{+})}\le {\lambda }_i\mu {\delta }_0{L}_{{Г}_0}^m{r}_1\underset{m\to \infty }{\to }0.\end{array}\right.$$

(2.29)

Therefore, based on the conclusions of Picard’s method, we obtain that the sequences ${\left\{{\theta }_m\right\}}_0^{\infty }$ and ${\left\{v_{i,m}\right\}}_0^{\infty }$ converge uniformly to the functions $\theta (x,t)$ and $v_i(x,t),(i=\overline{1,n})$ in $G_h^1(D_0),$ respectively, at that stage

$$\left\{\begin{array}{l}{\mathsf{\nu }}_i\in G_h^1(D_0):{‖{\mathsf{\nu }}_{\mathrm{i}}‖}_{G_h^1(D_0)}\le M_0,(0<M_0={\beta }_{11}+{\beta }_{12}=const;i=\overline{1,n}),\\ {\mathsf{\nu }}_{it}={\lambda }_i{\theta }_t+\exp (-{\displaystyle \frac{t}{\mu {\delta }_0}})\left\{\phantom{2}\right.-{\displaystyle \frac{1}{\mu {\delta }_0}}{J}_i+{\displaystyle \frac{1}{\sqrt{{\pi }^n}}}{\displaystyle \underset{R^n}{\overset{\phantom{2}}{\int }}(\exp }(-{\left|\xi \right|}^2))\sqrt{\mu }{\displaystyle \sum _{j=1}^n\frac{{\xi }_j}{\sqrt{t}}}{\psi }_{i{l}_j}(x+\\ +2\xi \sqrt{\mu t})d\left.\xi \right\},(l=x+2\xi \sqrt{\mu t}\in R^n),\\ {\displaystyle \sum _{0\le \left|k\right|\le 2}^{\phantom{2}}{‖D_{\phantom{2}}^k{\mathsf{\nu }}_{\mathrm{i}}‖}_{C(D_{+})}}\le {\beta }_{11};{‖{\mathsf{\nu }}_{it}‖}_{1h}\le {\beta }_{12},(i=\overline{1,n}).\end{array}\right.$$

(2.30)

Hence, based on condition (2.30) it follows

$${‖\mathsf{\nu }‖}_{G_{n,h}^1(D_0)}={\displaystyle \sum _{i=1}^n{‖{\mathsf{\nu }}_i‖}_{G_h^1(D_0)}=}{\displaystyle \sum _{i=1}^n\{}{\displaystyle \sum _{0\le \left|k\right|\le 2}^{\phantom{2}}{‖D^k{\mathsf{\nu }}_i‖}_{C(D_{+})}}+\underset{x\in R^n}{\sup }{\displaystyle \underset{0}{\overset{\infty }{\int }}h(t)\left|{\mathsf{\nu }}_{it}(x,t)\right|dt}\}\le N_0=const.$$

(2.31)

As a result, we obtain the following:

Theorem 1. Under the assumptions of Lemmas 1 and 2, and condition (2.28), the Cauchy problem (2.14) has a unique solution in $G_h^1(D_0).$ Furthermore, taking equation (2.2) into account, the n-dimensional Navier–Stokes system (1.1) with conditions (1.2) and (1.3) admits a unique solution in $G_{n,h}^1(D_0)$.

Note 2. Under the conditions of Theorem 1, the solution to the n-dimensional Navier–Stokes problem (1.1)–(1.3), obtained via rule (2.2), possesses second-order smoothness with respect to the spatial coordinates and first-order partial derivatives with respect to the time variable $t$, when $t>0$ (the conditionally smooth solution), since the function $\nu \in G_{n,h}^1(D_0)$.

Thus, condition (a.4) is proven.

3. Conclusions

The main contribution of this work is the development of a method that, without the need for additional assumptions, equivalently transforms the n-dimensional NS problem into a non-homogeneous linear equation of the heat conduction type with variable coefficients. In this context, the proposed algorithm (2.2) plays a crucial role in linearizing the n-dimensional NS system (1.1). First, by using formula (2.2), the nonlinear convective terms in the NSE are linearized with respect to the functions $\theta ,{\theta }_{x_i},(i=\overline{1,n})$. Second, by considering equation (2.2), we obtain Poisson-type equations for the pressure. Third, by eliminating the pressure from equation (2.9), we derive a linear parabolic problem (2.14), which is equivalently transformed into a regular second-kind boundary value problem with respect to the viscosity $\mu \in (0,1)$. Therefore, the unique solvability of the NSE problem for incompressible fluid with viscosity in $G_{n,h}^1(D_0)$ follows from the unique solvability of the resulting system of integral equations based on (2.2).

It is worth noting that even if the known functions ${{\rm Y}}_i,(i=1,2)$ are continuous but their partial derivatives are understood in the Sobolev sense [12], the algorithm (2.2) is still applicable. This fact is also one of the key advantages of the proposed method.