Reduced basis approximation and a posteriori error bounds for 4D-Var data assimilation

Abstract

We propose a certified reduced basis approach for the strong- and weak-constraint four-dimensional variational (4D-Var) data assimilation problem for a parametrized PDE model. While the standard strong-constraint 4D-Var approach uses the given observational data to estimate only the unknown initial condition of the model, the weak-constraint 4D-Var formulation additionally provides an estimate for the model error and thus can deal with imperfect models. Since the model error is a distributed function in both space and time, the 4D-Var formulation leads to a large-scale optimization problem for every given parameter instance of the PDE model. To solve the problem efficiently, various reduced order approaches have therefore been proposed in the recent past. Here, we employ the reduced basis method to generate reduced order approximations for the state, adjoint, initial condition, and model error. Our main contribution is the development of efficiently computable a posteriori upper bounds for the error of the reduced basis approximation with respect to the underlying high-dimensional 4D-Var problem. Numerical results are conducted to test the validity of our approach.

Appendix 1: Continuous 4D-Var formulation

The strong-constraint 4D-Var problem for a linear parabolic PDE on $(0,T)\times {\varOmega }$, with ${\varOmega }$ a Lipschitz domain,

$$\begin{aligned} \partial _ty + Ay = f, \quad y|_{\partial {\varOmega }}=0, \qquad y(0)=u, \end{aligned}$$

(68)

is classically rewritten as the optimal control problem:

$$\begin{aligned} \text {Find } (y^*,u^*) \in \mathop {\hbox {arginf}}_{ (y,u)\in \mathcal {Y}\times U \text { satisfies (68) } } J \end{aligned}$$

with a lower semi-continuous cost functional

$$\begin{aligned} J = \frac{\lambda }{2} \left||u-u_d\right||^2_U + \frac{1}{2} \int _0^T \left||Cy-z_d\right||^2_D \end{aligned}$$

(69)

on the tensor-product of $\mathcal {Y} := \{y\in L^2(0,T;H^1_0({\varOmega })); \partial _t y\in L^2(0,T;H^{-1}({\varOmega }))\}$ and $U := L^2({\varOmega })$. If the observation operator C has a unique continuation in $\mathcal {Y}$, J is coercive and strictly convex. Then, if $f\in L^2((0,T)\times {\varOmega })$ so the set of admissible states is non-empty, there exists a unique solution, see e.g. Fursikov (2000). To characterize and compute the solution, one can use duality techniques following Pontryagin et al. (1964) or Ekeland and Temam (1976). On introducing a Lagrange multiplier $(p^*,v^*)\in \mathcal {Y}\times U$, $p^*(T)= 0$ for the constraint, it is classical that the solution should satisfy (Marchuk and Shutyaev 2002)

$$(u^*,\psi )_U-(\lambda v^*,\psi )_U=(u_d,\psi )_U \ \forall \psi \in U$$

(70a)

$$(p^*(0),\varphi )_Y - (v^*,\varphi )_U= 0$$

(70b)

$$(Cy^*,C\varphi )_D-\left( \partial _tp^*-A^Tp^*,\varphi \right) _Y=(Cz_d,C\varphi )_D \ \forall \varphi \in \mathcal {Y}$$

(70c)

$$(y^*(0),\phi )_Y - (u^*,\phi )_U= 0$$

(70d)

$$\left( \partial _ty^* + Ay^*,\phi \right) _Y= (f,\phi )_Y \ \forall \phi \in \mathcal {Y}$$

(70e)

which is a well-posed saddle-point problem, well-approximated by the discretization (8) (Ern and Guermond 2010), again under the condition that the observation operator C has a unique continuation in $\mathcal {Y}$. Note that first adequately discretizing J then leads to exactly the same discrete Euler–Lagrange equations as (8).

The weak-constraint 4D-Var problem is also classical, see e.g. Fursikov (2000). The optimal control problem becomes

$$\begin{aligned} \text {Find } (y^*,u^*) \in \mathop {\mathrm{arginf}}_{ (y,u)\in \mathcal {Y}_{y_0}\times \mathcal {U} \text { satisfies (71) } } J \end{aligned}$$

for

$$\begin{aligned} \partial _ty + Ay = f + Bu, \quad y|_{\partial {\varOmega }}=0, \qquad y(0)=y_0, \end{aligned}$$

(71)

with the lower semi-continuous cost functional

$$\begin{aligned} J = \frac{1}{2} \int _0^T \left||u-u_d\right||^2_U + \frac{1}{2} \int _0^T \left||Cy-z_d\right||^2_D \end{aligned}$$

(72)

on the tensor-product of $\mathcal {Y}_{y_0} := \{y\in L^2(0,T;H^1_0({\varOmega })); \partial _t y\in L^2(0,T;H^{-1}({\varOmega })); y(0)=y_0 \}$ and $\mathcal {U}:= L^2((0,T)\times {\varOmega })$; the saddle-point becomes

$$\begin{aligned} (u^*,\psi )_U-(B^Tp^*,\psi )_Y&=(u_d,\psi )_U \ \forall \psi \in U \end{aligned}$$

(73a)

$$\begin{aligned} (Cy^*,C\varphi )_D-\left( \partial _tp^*-A^Tp^*,\varphi \right) _Y&=(Cz_d,C\varphi )_D \ \forall \varphi \in \mathcal {Y} \end{aligned}$$

(73b)

$$\begin{aligned} \left( \partial _ty^* + Ay^*,\phi \right) _Y - (Bu^*,\phi )_Y&= (f,\phi )_Y \ \forall \phi \in \mathcal {Y} \end{aligned}$$

(73c)

while existence and uniqueness of a soluion still hold under the same conditions.