Continuous data assimilation for the 3D primitive equations of the ocean

In this article, we show that the continuous data assimilation algorithm is valid for the 3D primitive equations of the ocean. Namely, the $L^2$ norm of the assimilated solution converge to that of the reference solution at an exponential rate in time. We also prove the global existence of strong solution to the assimilated system.


Introduction
In this paper, we address the continuous data assimilation algorithm applied to the 3D primitive equations of the ocean. We first show that the data assimilation primitive system possesses a global in time strong solution. Then, we prove that the strong solution of the assimilated system converges to that of the reference system exponentially fast in both L 2 norm. Namely, the continuous data assimilation algorithm, which we elaborate in details below, is valid for the 3D primitive equations of the ocean. This work sets up the analytical foundation and provides the guideline for the application of the new and promising continuous data assimilation algorithm for the primitive equations of the ocean, which are widely used in geophysics and meteorology.
We start with an introduction to the idea of data assimilation following [67]. Simulations for dynamical systems, such as the primitive equations, require accurate initialization process in order to make reliable predictions. The initialization procedure depends on how acquired observations, such as measurements of temperature and velocity, are properly interpolated in space and time, to complete the information throughout the entire spacetime domain, while maintaining the dynamical balance between these fields. Even if direct observations of these fields are available, typically they are not uniformly distributed in space or are very sparse. Moreover, the errors contained in the measurements and parameters of the model, combined with the highly nonlinear dynamics of the PDEs, makes basic interpolation unfitting for initial condition construction or selection, even for short term predictions. Meteorologists have designed a series of diagnostic experiments when creating precise initializing processes that minimize the loss of information from the measurements. For instance, the interpolated function in space and time must satisfy the conservation laws in the continuous model Date: May 17, 2018. equations. Meteorologists also use various combinations of collected information about the state of the system to see which combination(s) will drive the system closer to the collected data. In this process, they may also use other forms of collected measurements such as, the average temperature of a small region obtained through image processing of satellite observation. One must know how to properly make use of this information, in particular when the collected data is not one of the evolving state variable in the dynamical system. In the context of meteorology and atmospheric physics, data assimilation algorithm when some state variable observations are not available as an input, has been studied in [1][2][3]67] and the references therein, for some simplified numerical forecast models. Although nonlinear interactions between the scales of motion and the parameters in the system exist, it has been shown in several settings that if the dynamical model used as a source of a priori information captures the important properties of the system being modeled, then, one can identify the full state of the system knowing only coarse observation from partial data that is selected properly.
In this work, we consider the continuous data assimilation introduced by Azouani, Olson, and Titi in [4,5] (see also [6][7][8] for early ideas in this direction). This approach, which we call AOT data assimilation or the linear AOT algorithm, is based on feedback control at the partial differential equation level. We now describe the idea of the AOT algorithm. Consider a dynamical system in a general form, For example, this could represent a system of partial differential equations modeling fluid flow in the atmosphere or the ocean. A central difficulty is that, even if one were able to solve the system exactly, the initial data u 0 is largely unknown. For example, in a weather or climate simulation, the initial data may be measured at certain locations by weather stations, but the data at locations in between these stations remain unknown. Therefore, one might not have access to the complete initial data u 0 , but only to the observational measurements, which we denote by I δ (u 0 ). (Here, I δ is assumed to be a linear operator that can be taken, for example, to be an interpolation operator between grid points of maximal spacing h, or as an orthogonal projection onto Fourier modes no larger than k ∼ 1/h.) Moreover, the data from measurements may be streaming in moment by moment, so in fact, one often has the information I δ (u) = I δ (u(t)), for a range of times t. Data assimilation is an approach that eliminates the need for complete initial data and also incorporates incoming data into simulations. Classical approaches to data assimilation are typically based on the Kalman filter. See, e.g., [61][62][63] and the references therein for more information about the Kalman filter.
In 2014, the AOT algorithm was introduced in [4,5]. This new approach overcomes some of the drawbacks of the Kalman filter approach (see, e.g., [9] for further discussion). The approach has been the subject of much recent study in various contexts, see, e.g., [10-21, 68, 69]. Specifically, the following system was proposed and studied in [4,5]: This system, used in conjunction with (1.1), forms the AOT algorithm for data assimilation of system (1.1). In the case where the dynamical system (1.1) is the 2D Navier-Stokes equations, it was proven in [4,5] that, for any divergence-free initial dataũ 0 ∈ L 2 , u(t) −ũ(t) L 2 → 0, exponentially in time. In particular, even without knowing the initial data u 0 , the solution u can be approximately reconstructed for large times. We emphasize that, as noted in [4], the initial data for (1.2) can be any L 2 function, even v 0 = 0. Thus, no information about the initial data is required to reconstruct the solution asymptotically in time.
In the description below we denote the two-dimensional horizontal gradient, Laplacian, and divergence operators by ∇ 2 , ∆ 2 , and div 2 , respectively. The 3D primitive equations of the ocean and atmosphere (PEs) on the domain Ω = R 2 /Z 2 × (−h, 0) with h > 0 read: where v = (v 1 , v 2 ) and w stand for the horizontal and vertical components of the three-dimensional velocity field u = (v, w) and θ for the temperature and p for the scalar pressure, all of which are unknowns of the equations. Also, f is the Coriolos parameter and q is the given heating source. The viscosity and the heat diffusion operators L 1 and L 2 are given by where the positive constants ν v , ν w , and η v , η w are horizontal and vertical viscosities and heat diffusivities, respectively. For simplicity, we assume that the initial data is sufficiently smooth (made more precise below) and satisfy The primitive equations of the ocean are considered to be the fundamental model for meteorology and climate prediction (c.f. [64]). The primitive equations are derived from the full compressible Navier-Stokes equations, which are too complicated in terms of application and contain phenomena such as shocks and sound waves that are not interesting from the perspective of geophysics. Thus, by applying the Boussinesq and the hydrostatic approximations to the full Navier-Stokes system, one obtains the primitive system consists of the horizontal momentum equations, the conservation of mass, and the thermodynamic equation-diffusion of temperature.
From the point of view of data assimilation, the primitive equations are ideal in applying the algorithm since the solution to the system has chaotic large-time behavior; and it also has a finite-dimensional global attractor (see [54][55][56] and the reference therein) making PE an excellent candidate for studying large-time behavior. Therefore, in view of the AOT data assimilation, we have: where 1/δ is the spatial resolution and the operator I δ satisfies the following inequality: In particular, one may take I δ to be orthogonal projection onto the first c/δ Fourier modes, for some constant c. Remark 1.1. Other physically relevant choices of I δ , such as a nodal interpolation operator, have been considered in the case of the linear AOT algorithm (see, e.g., [19]).
The principal aim of this article is double-fold. First, we show the existence of global strong solution to the assimilation system, which is necessary for the convergence results. For the proof, we follow the idea in the classic paper [30] and provide the a priori estimates for the L 2 , L 6 , and H 1 norms of the solutions. The main difficult comes from the lack of the momentum equation for the third component of the velocity. We use the anisotropic estimates as well as the divergence-free condition to bound the nonlinear terms. Secondly, to the best knowledge of the author, we present the first continuous data assimilation result for the full primitive equations of the ocean in 3D. We prove that with a sufficiently refined resolution (δ small) and large enough nudging coefficients (µ u and µ θ ), the assimilated solution approaches the reference solution at an exponential rate in the sense of L 2 norm.

Preliminaries and Summary of Main Results
In this paper, we frequently use the L 2 norm and H 1 semi-norms, defined by All through this paper we denote ∂ j = ∂/∂x j , ∂ jj = ∂ 2 /∂x 2 j , ∂ t = ∂/∂t, and the horizontal gradient ∇ 2 = (∂ 1 , ∂ 2 ) and horizontal Laplacian ∆ 2 = ∂ 11 + ∂ 22 . Also, we denote the usual Lebesgue and Sobolev spaces by L p x and H s x ≡ W s,2 x , respectively, and let F be the set of all trigonometric polynomial over T 3 and define the space of zero-average smooth functions Also, we define H and V to be the closures of V in L 2 x and H 1 x , with inner products respectively, associated with the norms |u| = (u, u) 1/2 and u = ((u, u)) 1/2 . The latter is indeed a norm due to Poincaré inequality and denote the closure of V θ in L 2 and H 1 by H θ and V θ , for θ. We refer to readers to [65,66] for more standard details on Navier-Stokes, which will be used in this paper.
The following Gagliardo-Nirenberg-Sobolev type inequality is frequently used in our estimates.
We list some special cases below, and the constant C may vary but is independent of φ and u, respectively.
The next lemma about anisotropic estimates is critical to out treatment of the nonlinear terms.
Now we state our main results of this paper. The first theorem concerns the global existence of solution to system (1.4).
The next theorem provides the exponential-in-time convergence in L 2norm of the continuous data assimilation for the 3D primitive equations of the ocean.
where C depends on λ 1 , ν v , ν w , η v , η w , the domain Ω, and the reference solution, and with Recall the definition of the strong solution to (1.3) as follows, In order to prove our main theorems, we need the following theorem summarized from the results in [30] on the global well-posedness of system (1.3).
Theorem 2.7. Let Q ∈ H 1 , v 0 ∈ V , θ 0 ∈ V θ and T > 0. Then, there exists a unique strong solution of (1.3) on [0, T ] which continuous depends on the initial data.
Further, we summarize the following theorems from [54][55][56] regarding the uniform bounds of the L 2 , H 1 , and H 2 norms of the solution and the existence of finite dimensional global attractors of system (1.3). See also [70] for a new proof.
Moreover, there exists a bounded absorbing ball for all solutions (v, θ) of (1.3) in V × V θ . Theorem 2.9. With the same assumption as in Theorem 2.8, for the strong For the sake of completeness, we state a special case of the generalized Grönwall's inequality used in [17] and [57].
Lemma 2.10. Suppose that Y (t) is a locally integrable and absolutely continuous function that satisfies the following: Then, Y (0) → 0 at an exponential rate as t → ∞.

Existence of solution to the assimilated system
In this section, we provide the proof for the existence of solutions to system (1.4). We provide the a priori estimates only and one can obtain the solution via standard Galerkin approximation henceforth, which we omit for simplicity.
3.1. L 2 -estimates. First we multiplyṽ andθ to the corresponding equations in system (1.4), respectively, integrate over Ω, and get and without loss of generality, we assumep 0 = 0. Due to the divergence-free condition, the first two terms on the right side of the first equation above add up to zero. The third term is bounded by while for the assimilation term, we observe that where we used condition (1.5) and Poincaré's inequality. Thus, the last term is bounded by Now choose δ small enough so that where we used Cauchy-Schwarz inequality. Regarding the estimates of the right side ofθ-equation, we also observe that divergence-free condition implies that the first two terms sum up to zero, while the third term is bounded by Ω |q||θ| ≤ q L 2 θ L 2 .
As for the last term, similar to the case inṽ, we bound it by By choosing small enough δ such that where we used the fact that θ 2 L 2 ≤ 2h 2 ∂ zθ 2 L 2 + 2h θ (z = 0) 2 L 2 . Therefore, by Theorem 2.8 and Grönwall's inequality, we conclude that where K 1 depends only on θ H 1 and q L 2 . Whence, we obtain ṽ L 2 ≤ K 2 , and here K 2 depends on K 1 and v H 1 . Furthermore, we have the bounds for any T > 0.
3.2. L 6 -estimates. We begin by multiplying |ṽ| 4ṽ and |θ| 4θ to the equations ofṽ andθ, respectively, integrate over Ω, and get where we used the fact that f k ×ṽ|ṽ| 4 = 0. Thanks to the divergence-free condition ∂ z W = −∇ 2 ·ṽ, the first and second terms on the right side of both equations add up to zero, respectively. Then, we estimate the remaining terms on the right side of the above two equations. For the third term on the right side ofṽ-equation, we use Hölder's inequality and obtain The fourth term is bounded by where we used condition (1.6) for the operator I δ and Cauchy-Schwarz inequality. Next, we estimate the third and fourth terms on the right side of theθ-equation. By Hölder's inequality we first bound the third term by while regarding the last term, we follow similarly the estimates forṽ and bound it by Thanks to Grönwall's inequality, by summing up the above estimates and by Theorem 2.7 and Theorem 2.8, we get where K 3 depends on θ H 1 , q H 1 , and c 1 , δ, µ θ , λ 1 , while K 4 depends on v H 1 and K 3 .

H 1 -estimates.
In order to obtain the H 1 -bounds onṽ, we first need to estimate ∂ zṽ L 2 and ∇∂ zṽ L 2 .
3.3.1. ∂ zṽ L 2 estimates. Multiply the first equation in (1.4) by −∂ 2 zzṽ , integrate over Ω, and we get where we used the equation ∂ zp = −θ and the fact that Then, we estimate the four terms on the right side of the above equation. First, by Hölder's inequality and Lemma 2.3, we have where we integrated by parts in the first equality. Next, in order to bound the second term, we use the divergence-free condition and proceed as Thus, we get For the third term, integration by parts yields where we used Hölder's and Cauchy-Schwarz inequalities. Regarding the last term, by condition (1.6) of I δ and Poincaré's inequality, we have Hence, by combining all the above estimates and using Grönwall's inequality, we obtain where K 5 depends on δ, µ u and K 1 through K 4 .
We bound the first term on the right side of the above equation by Using divergence-free condition, we estimate the second term as where we used anisotropic Sobolev inequality. By Cauchy-Schwarz inequality, we bound the third term by As for the last term, by condition (1.6), we have where we used Hölder's and Minkowski inequalities. Combining all the above estimates, and by Grönwall's inequality, we obtain ∇ 2ṽ (T ) 2 L 2 < K 6 as well as for arbitrary T > 0, where K 6 depends on δ, µ u and K 1 through K 5 ..

θ H 1 estimates.
Multiplying theθ-equation in (1.4) by −(∆ 2 +∂ 2 zz )θ, integrate over Ω, we get where we used ∂ zw = −∇ 2 ·ṽ. Using Lemma 2.3, the first term is bounded by In order to estimate the second term, we apply the anisotropic estimates, i.e., Sobolev inequality in two-dimension, and bound it by By Cauchy-Schwarz inequality, the third term on the right side of the above equation is bounded by Regarding the last term, we use condition (1.6) and Hölder's and Minkowski inequalities and get Summing up all the above estimates, and using the H 1 estimates forṽ, as well as Grönwall's inequality, we get for arbitrary T > 0, where K 7 depends on q H 1 and K 1 through K 6 . Now, with the above L 2 , L 6 , and H 1 estimates ready, suppose that there exists some time T * > 0 such that the interval (0, T * ) is the maximal interval of existence of the solution (ũ,θ,p). We infer that must be valid provided T * < ∞. Therefore, we obtain a contradiction to the above estimates and the solution thus exists globally for all T > 0. Proof of Theorem 2.4 is complete.
4. Exponential-rate Convergence to the Reference Solution in L 2 -norm In this section, we prove Theorem 2.5.
Proof. By subtracting the corresponding equations of system (1.4) from those of system (1.3), we obtain where we denote by U = (V, W ) the difference V = v −ṽ, W = w −w, Θ = θ −θ, and P = p −p. For −h < z < 0, integrating the equation of the pressure P vertically, we get and without loss of generality, we assume P 0 = 0. Then, we multiply the equations of V and Θ by V and Θ, respectively, use the above equation of P , integrate over Ω, and get where we used the divergence-free condition ∂ z W = −∇ 2 · V . Next, we estimate the right sides of the above two equation of V and Θ. Once again in view of the divergence-free condition, we deduce that the third and fourth terms on the right side of the V -equation add up to zero. For the same reason, the third and fourth terms on the right side of the Θ-equation also sum up to 0. In order to bound the first term on the right side of the V -equation, we use Hölder's inequality and Lemma 2.3, and obtain For the second term on the right side of the V -equation, we estimate after integration by parts as Regarding the last but second term on the right side of the V -equation, we have where we integrated by parts and used Cauchy-Schwarz inequality. The last term in the V -equation is bounded as Next, we estimate the right side of the Θ-equation. For the first term, by Hölder's inequality and Lemma 2.3, we get In order to estimate the second term, we first integrate by parts, then proceed as As for the last term, we bound it as Combining all the above estimates, and by denoting and Y (t) = ∇V 2 L 2 + ∇Θ 2 L 2 , and using Theorem 2.7 and Theorem 2.8, we obtain d dt X(t) + 1 2 Y (t) (4.1) ≤ C + C ν,η ∇ 2 v 4 L 2 + ∆v 4 L 2 + ∇θ 4 L 2 + ∆θ 2 L 2 + ∆θ 4 where µ min = min{µ u , µ θ } and C depends on the reference solution, λ 1 , and the domain. Now choose where µ max = max{µ u , µ θ }, and choose such that the coefficient of X(t) on the right side of (??) is negative, where the constant C ν,η also depends on ν v , ν w , η v , η w . Then, by Grönwall's inequality, we conclude that X(t) tends to zero exponentially fast as time t → ∞. Thus, the exponential L 2 -convergence of the data assimilation solution to the reference solution is proven.
Remark 4.1. In view of the available results regarding the regularity of solutions to the 3D primitive equations and the continuous data assimilation on other dynamical systems, we believe that the H 1 norms of V and Θ should also tend to zero at an exponential rate and the relevant work will be presented in future works.